--- jupytext: text_representation: extension: .md format_name: myst kernelspec: display_name: Python 3 language: python name: python3 --- ```{contents} :local: :depth: 2 ``` # Problems: High-dimensional Classification :::{admonition} Get this problem set :class: seealso {download}`Download everything (Topic3.4-High-dimensional_Classification_Problems.zip) ` — the notebook and `membranes_fabrication.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. ::: :::{admonition} Before you start :class: tip This problem set accompanies {doc}`High-dimensional Classification `. 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 `data/membranes_fabrication.csv` is a gas-separation membrane database: 1,672 records covering more than 280 membranes, hand-extracted from 42 papers published between 2009 and 2022 ([github.com/tenoriolms/database_membranes](https://github.com/tenoriolms/database_membranes), Apache-2.0). Each record says what the membrane was made of, how it was mounted, and the conditions it was tested at. The task here is a **database audit** rather than a prediction anyone needs to make. Every record already states its material family. The question is whether a classifier can **recover that family from the fabrication and testing metadata alone**. If it can, the records are not independent samples of chemistry: they carry the fingerprint of how each sub-field runs and reports its experiments, and any model trained on the database inherits it. This is worth asking of a curated literature dataset before trusting anything built on it. The file has been prepared for that question. **No permeability or permeance column is in it.** The columns that would give the answer away outright are also gone: `subtype` (47 values, not one shared between two families, so it is the label spelled differently), `filler_loading` (non-zero for every mixed-matrix membrane and for nothing else), and the material description. ### `source_ref`, and why every split is by paper `source_ref` names the paper each row came from. It is **not a feature** — it is the **grouping key**, and every split below depends on it. The 1,672 rows came from 42 papers, and one paper contributes 134 of them: its own handful of membranes, measured at many conditions. Split those rows at random and rows from the same paper land on both sides, so a classifier can score well by recognizing the paper rather than the chemistry. Every split below is therefore made **by paper**, with `GroupShuffleSplit`, and no paper appears in both the training and the test set. Part C measures what that costs. Two more preparation decisions, both forced by the data: - **Four families, not six.** `CMS` (21 rows) and `ceramic` (24 rows) each come from a *single* paper, so a split by paper puts all of either class entirely on one side and the class becomes untrainable. They are dropped, leaving polymer, zeolite, MMM and MOF. - **Missing values are imputed with the training-set median**, inside the pipeline. The numeric columns are 16 to 99% complete — literature compilations record whatever each paper happened to report — and `SVC` cannot take a `NaN`. Doing it inside the pipeline is what keeps the test rows out of the median. ```{code-cell} ipython3 %matplotlib inline import numpy as np import pandas as pd import matplotlib.pyplot as plt from sklearn.model_selection import GroupShuffleSplit, GroupKFold, GridSearchCV from sklearn.pipeline import Pipeline from sklearn.compose import ColumnTransformer from sklearn.preprocessing import OneHotEncoder, StandardScaler from sklearn.impute import SimpleImputer from sklearn.svm import SVC from sklearn.metrics import confusion_matrix, classification_report try: plt.style.use('../settings/plot_style.mplstyle') # available inside the book except OSError: pass # downloaded notebook: use defaults df_all = pd.read_csv('data/membranes_fabrication.csv') FAMILIES = ['polymer', 'zeolite', 'MMM', 'MOF'] df = df_all[df_all['material_type'].isin(FAMILIES)].reset_index(drop=True) NUMERIC = ['temperature_K', 'feed_pressure_Pa', 'permeate_pressure_Pa', 'delta_pressure_Pa', 'feed_flow_rate_mL_min', 'sweep_gas_flow_mL_min', 'stage_cut_pct', 'min_thickness', 'max_thickness', 'min_pore_size_nm', 'max_pore_size_nm', 'total_pore_volume_cc_g', 'micropore_volume_cc_g', 'specific_surface_area_m2_g', 'surface_area_cm2', 'aging_days'] CATEGORICAL = ['configuration', 'support_material', 'pore_size_type', 'sweep_gas', 'provided_data_type'] X = df[NUMERIC + CATEGORICAL] y = df['material_type'].values groups = df['source_ref'].values train_idx, test_idx = next( GroupShuffleSplit(n_splits=1, test_size=0.3, random_state=0).split(X, y, groups=groups)) X_train, X_test = X.iloc[train_idx], X.iloc[test_idx] y_train, y_test = y[train_idx], y[test_idx] groups_train = groups[train_idx] print(f"{len(df_all)} records, {len(df)} in the four retained families") print(df['material_type'].value_counts().to_string()) print(f"\n{len(train_idx)} training rows from {len(set(groups_train))} papers") print(f"{len(test_idx)} test rows from {len(set(groups[test_idx]))} papers") print(f"no paper on both sides: {set(groups_train).isdisjoint(set(groups[test_idx]))}") print(f"test composition: {pd.Series(y_test).value_counts().to_dict()}") ``` One helper builds the pipeline, with the categorical branch switchable. Parts A and B both turn on that switch. ```{code-cell} ipython3 def make_pipeline(use_categorical, C=10.0, gamma=0.1, categorical=None): """RBF-kernel SVC on the numeric columns, optionally plus one-hot categoricals. The numeric branch imputes with the training median and then standardizes; both are fitted inside the pipeline, so the test rows never contribute to either. """ cats = CATEGORICAL if categorical is None else categorical branches = [('num', Pipeline([('impute', SimpleImputer(strategy='median')), ('scale', StandardScaler())]), NUMERIC)] if use_categorical and cats: branches.append(('cat', OneHotEncoder(handle_unknown='ignore', sparse_output=False), cats)) return Pipeline([('prep', ColumnTransformer(branches)), ('clf', SVC(kernel='rbf', C=C, gamma=gamma))]) def n_features(pipe): """How many columns the fitted preprocessor actually produces.""" return pipe.named_steps['prep'].transform(X_train).shape[1] print("make_pipeline ready") ``` `handle_unknown='ignore'` is not optional here. The test papers are papers the encoder has never seen, so they can easily report a support material or a sweep gas that appears nowhere in the training rows; without it, `transform` would raise. 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. ```{code-cell} ipython3 :tags: [skip-execution] # Run this once. It defines `check`, used after each Part A answer below. import hashlib _EXPECTED = { 'q1': (4, ['09cb69a6e0d5', '81d0971f19d3', 'd72c71286bf8', '192160bf535f', 'ebc96010c6fc', '02d44139d45b', '167cd7e4f245', 'b7ba25a2f627', '8bb947280191', '8cc2062b35c6', '93473885b84e']), 'q2': (4, ['4f62466d8ec8', 'c40b0be84707', '88c288ebcc67', '06e0abe23cef', '143466381d8e', '359b7076a95a', 'eaee3e6ae664', '0f121a5d3e10', '800592cbbeef', '43f765f03bf7', '76bb34025409', '953800322c3e', '6e4a5cbbb332']), 'q3': (4, ['4b16dee212a6', 'd1ccb042876a', '0d58cbf3cc90', 'f0613dc0d593', '5d5d32037dc9', '87eb2ec0799b', '0895593323f9', '6d400c5ae1e3', '1a5b887be48d', 'bd1e16883bbc', '95bb49df2c33', 'd72aadfd7bcd', 'cd61eaf3e16f', '866fb79d1805', '8283afb503f3']), } 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'}") ``` --- ## Part A — Skill Checks (30 pts) Three questions, 10 points each, and together they are one controlled experiment. A1 is the number a classifier has to beat. A2 and A3 then use **the same split, the same kernel and the same hyperparameters**, and differ only in whether the five categorical columns are one-hot encoded and handed to the model. The size of that difference is what Parts B and C are about. Leave `C=10` and `gamma=0.1` alone in this part — tuning them is Part C's job. Fit on the training rows, score on the test rows. ### A1. The number to beat (10 pts) :::{exercise} :label: pr-cls-mem-baseline A classifier that ignores every column and always answers with the most common family still scores something, and every model below has to beat it. Using `y_test` only, compute the fraction of test rows belonging to the **largest** of the four families. Assign it to `baseline_accuracy`. ::: ```{code-cell} ipython3 :tags: [skip-execution] # YOUR CODE HERE baseline_accuracy = ... ``` ```{code-cell} ipython3 :tags: [skip-execution] check("q1", baseline_accuracy) ``` ### A2. Numeric columns only (10 pts) :::{exercise} :label: pr-cls-mem-numeric Build the pipeline with `use_categorical=False`, so the model sees only the 16 numeric fabrication and operating columns. Fit it on the training rows and score it on the test rows. Assign the test accuracy to `acc_numeric`. ::: ```{code-cell} ipython3 :tags: [skip-execution] # YOUR CODE HERE acc_numeric = ... ``` ```{code-cell} ipython3 :tags: [skip-execution] check("q2", acc_numeric) ``` ### A3. Adding the categorical columns (10 pts) :::{exercise} :label: pr-cls-mem-onehot Now build the same pipeline with `use_categorical=True`, so the five categorical columns are one-hot encoded and appended to the standardized numeric ones. Nothing else changes — same split, same kernel, same `C` and `gamma`. Fit, score on the test rows, and assign the test accuracy to `acc_onehot`. Print the feature count both ways while you are here. Five columns became 20. ::: ```{code-cell} ipython3 :tags: [skip-execution] # YOUR CODE HERE acc_onehot = ... ``` ```{code-cell} ipython3 :tags: [skip-execution] check("q3", acc_onehot) ``` --- ## Part A summary | | Features | Test accuracy | |---|---|---| | Majority family | 0 | A1 | | Numeric only | 16 | A2 | | Numeric + one-hot | 36 | A3 | --- ## Part B — Visualization (35 pts) In Part A, one-hot encoding gained more than every numeric column put together. That result is unspecific: five columns went in at once, so nothing yet says *which* of them carries the information, or which families it separates. This part pins the effect down in two directions. Per-class metrics say **which families** the encoding rescued — overall accuracy hides that one family was already being classified well while another was barely being found at all. An ablation says **which column** does the work. Keep two groups of columns apart as you write: the **material** categoricals describe what the membrane is made of and how it is mounted (`support_material`, `configuration`), while the **protocol** categoricals describe how the permeation test was run and reported (`sweep_gas`, `pore_size_type`, `provided_data_type`). Which group wins is the thing to report. :::{exercise} :label: pr-cls-mem-ablation **1. Compare the two models class by class.** Produce a grouped bar chart of **per-class recall** for the numeric-only and the one-hot model side by side, four families on the $x$-axis. `classification_report(..., output_dict=True)` or `precision_recall_fscore_support` will give you the numbers. Print the full classification report for both models as well. **2. Show both confusion matrices.** One figure, two panels, same family order on both axes, counts annotated in the cells. Use the same color scale in both panels so they can be read against each other. **3. Ablate the categorical columns.** For each of the five, fit the pipeline with **that column as the only categorical** (`make_pipeline(True, categorical=[col])`) and record the test accuracy. Plot the five as a bar chart against the numeric-only accuracy as a reference line. **4. Look at the winner's contents.** Take the categorical column that gained the most in task 3 and cross-tabulate it against `material_type` over all 1,627 rows (`pd.crosstab`). Show it as an annotated heatmap of **row-normalized** shares, so each row reads as $P(\text{family} \mid \text{level})$. Then answer, in **one or two sentences each**: - **(a)** Which family gained the most recall from one-hot encoding, and which family was already being classified well without it? - **(b)** Which single categorical column gained the most, and is it a *material* column or a *protocol* column? Name which one and what that implies about what the model is reading. - **(c)** Pick one row of your task-4 heatmap that is almost entirely one family, and say in physical terms why that experimental choice goes with that material family. Label your axes. Annotate the confusion-matrix cells with counts. ::: ```{code-cell} ipython3 :tags: [skip-execution] # YOUR CODE HERE ``` --- ## Part C — Open Ended (35 pts) Everything so far has used `C=10` and `gamma=0.1`, values handed to you with no justification, and a split by paper presented as simply the right thing to do. Part C tests both choices. The second one turns out to matter far more than the first. The chapter's opening section gives the correct procedure: split first, tune inside the training set, touch the test set once. That protocol has a hole in it when rows are not independent. You can follow every step and still leak, because the leak is in the split itself rather than in the tuning. With 1,672 rows drawn from 42 papers, this database is exactly that case. Measure the effect, then say what the database can honestly support. :::{exercise} :label: pr-cls-mem-audit **1. Tune the kernel properly.** Run `GridSearchCV` over $C \in \{0.1, 1, 10, 100, 1000\}$ and $\gamma \in \{10^{-4}, 10^{-3}, 10^{-2}, 10^{-1}, 1\}$ for the one-hot pipeline, using `GroupKFold(n_splits=3)` on the **training rows only** and passing `groups=groups_train` to `fit`. Report the best $C$ and $\gamma$, the cross-validated score, and the test accuracy. Compare against the untuned A3 result and say how much tuning was worth. **2. Measure the leak.** Refit the same one-hot pipeline under a **random row-wise** split (`train_test_split(..., test_size=0.3, random_state=0, stratify=y)`) and report its test accuracy next to the grouped one. Report each split's own majority-class baseline too — they differ, so the raw accuracies are not directly comparable without them. **3. Test whether the leaky-looking column earns its place.** `support_material` looks like it should be a label in disguise: several of its levels contain exactly one family. Refit the full one-hot pipeline **without** it and report the test accuracy. The result is not what the crosstab predicts — explain the discrepancy in one or two sentences. **4. Say what the database supports.** In **two or three sentences**: state whether this database can be used to train a model that predicts membrane behavior for a *new* research group's membrane, and cite the two numbers from task 2 that decide it. There is more than one defensible position in task 4. You are graded on whether the numbers you produced support the one you take. ::: ```{code-cell} ipython3 :tags: [skip-execution] # YOUR CODE HERE ``` --- ## Summary - **Part A** ran one controlled experiment: the majority-family baseline, an RBF-kernel SVC on 16 numeric columns, and the same model with five categorical columns one-hot encoded into 20 more. The encoding was worth more than every numeric column combined. - **Part B** localized that gain — which families it rescued, and which single column carried it. The answer was a column describing the measurement protocol, not the material. - **Part C** tuned the kernel, found it barely moved the result, and then measured what the split by paper had been protecting against all along. ## Additional Reading 1. The database: [github.com/tenoriolms/database_membranes](https://github.com/tenoriolms/database_membranes) (Apache-2.0), and L. M. S. Tenorio et al., "Construction and analysis of a detailed new database for gas separation membranes", *Observatório de la Economía Latinoamericana* (2024). 2. scikit-learn User Guide: [Cross-validation for grouped data](https://scikit-learn.org/stable/modules/cross_validation.html#cross-validation-iterators-for-grouped-data), [SVM](https://scikit-learn.org/stable/modules/svm.html), [GridSearchCV](https://scikit-learn.org/stable/modules/grid_search.html), [Encoding categorical features](https://scikit-learn.org/stable/modules/preprocessing.html#encoding-categorical-features). 3. L. M. Robeson, "The upper bound revisited", *J. Membr. Sci.* **320**, 390 (2008). [doi:10.1016/j.memsci.2008.04.030](https://doi.org/10.1016/j.memsci.2008.04.030) — the permeability/selectivity trade-off that organizes this field, and the reason the permeability columns were withheld from this file.