""" PatientLeakageClassifier — runs the RF→GridSearchCV→SVM pipeline with patient-aware or image-level train/test splitting. """ import os, csv, re import numpy as np from sklearn.model_selection import ( train_test_split, GridSearchCV, StratifiedGroupKFold ) from sklearn.ensemble import RandomForestClassifier from sklearn.svm import SVC from sklearn.pipeline import Pipeline from sklearn.preprocessing import StandardScaler from sklearn.metrics import accuracy_score # --------------------------------------------------------------------------- # Helpers # --------------------------------------------------------------------------- def _f2n(fname): m = re.search(r'\((\d+)\)', fname) num = int(m.group(1)) if m else None for cls_key, prefix in [("Bengin cases","B"),("Malignant cases","M"),("Normal cases","N")]: if fname.startswith(cls_key.rstrip("s")): return f"{prefix}_{num:03d}" if num else fname return fname DEFAULT_NFEATURES = [50, 100, 200, 300, 400, 500, 750, 1000, 1500, 2000, 3000, 5000] DEFAULT_GAMMA_RANGE = (-1.5, 1, 7) # logspace args class PatientLeakageClassifier: """Run the classification pipeline for one model, seed, and split type. Parameters ---------- manifest_path : str Path to patient manifest CSV. features_dir : str Directory with {Model}_features.npz files. n_jobs : int Parallelism for RF and GridSearchCV. """ def __init__(self, manifest_path, features_dir, n_jobs=8): self.features_dir = features_dir self.n_jobs = n_jobs self._image_to_patient = {} with open(manifest_path, newline="") as f: reader = csv.DictReader(f) img_col = ("confirmed_images" if "confirmed_images" in reader.fieldnames else "images") for row in reader: pid = row["patient_id"] for img in row[img_col].split(";"): if img: self._image_to_patient[img] = pid # ------------------------------------------------------------------ # Public API # ------------------------------------------------------------------ def run(self, model_name, seed, split_type="patient", nfeatures_list=None, gamma_logspace=None, return_predictions=False): """Run full pipeline: train/test split → RF importance → GridSearchCV over nfeatures×gamma → retrain → test. Parameters ---------- model_name : str e.g. "VGG16" seed : int Random seed for the train/test split. split_type : str "image" or "patient". nfeatures_list : list Feature counts to search over. gamma_logspace : tuple Args for np.logspace (start, stop, n). return_predictions : bool If True, also include the test labels/predictions ("y_true", "y_pred") so callers can build a confusion matrix without refitting. Returns ------- dict with keys: model, seed, split_type, cv, test, nfeat, gamma, n_train, n_test, n_train_patients, n_test_patients """ if nfeatures_list is None: nfeatures_list = DEFAULT_NFEATURES if gamma_logspace is None: gamma_logspace = DEFAULT_GAMMA_RANGE X, Y, groups = self._load_data(model_name) (X_tr, X_te, y_tr, y_te, g_tr, cv, n_tr_patients, n_te_patients) = self._split(X, Y, groups, seed, split_type) sorted_idx = self._rf_ranking(X_tr, y_tr) best_cv, best_n, best_g = self._grid_search( X_tr, y_tr, g_tr, cv, sorted_idx, nfeatures_list, gamma_logspace) # --- Retrain + test --- final = Pipeline([("scaler", StandardScaler()), ("svm", SVC(kernel="rbf", C=10, gamma=best_g))]) final.fit(X_tr[:, sorted_idx[:best_n]], y_tr) y_pred = final.predict(X_te[:, sorted_idx[:best_n]]) result = { "model": model_name, "seed": seed, "split_type": split_type, "cv": float(best_cv), "test": float(accuracy_score(y_te, y_pred)), "nfeat": best_n, "gamma": float(best_g), "n_train": X_tr.shape[0], "n_test": X_te.shape[0], "n_train_patients": n_tr_patients, "n_test_patients": n_te_patients, } if return_predictions: result["y_true"] = [str(v) for v in y_te] result["y_pred"] = [str(v) for v in y_pred] return result def cv_curve(self, model_name, seed, split_type="image", nfeatures_list=None, gamma_logspace=None): """CV accuracy as a function of the number of RF-selected features. Same split → RF ranking → per-nfeatures γ-grid as ``run``, but returns the whole curve (used by the Figure 3 / S1 / S4 CV-vs-features plots) instead of only the best point. Returns ------- dict: {"curves": {nfeat: cv_acc}, "best_nfeat", "best_cv", "best_gamma"} """ if nfeatures_list is None: nfeatures_list = DEFAULT_NFEATURES if gamma_logspace is None: gamma_logspace = DEFAULT_GAMMA_RANGE X, Y, groups = self._load_data(model_name) X_tr, _, y_tr, _, g_tr, cv, _, _ = self._split( X, Y, groups, seed, split_type) sorted_idx = self._rf_ranking(X_tr, y_tr) pipe = Pipeline([("scaler", StandardScaler()), ("svm", SVC(kernel="rbf", C=10))]) valid_nfeats = [n for n in nfeatures_list if n <= X_tr.shape[1]] curves, best_cv, best_n, best_g = {}, 0, None, None for nfeat in valid_nfeats: pg = {"svm__gamma": (1.0 / nfeat) * np.logspace(*gamma_logspace)} grid = GridSearchCV(pipe, pg, cv=cv, scoring="accuracy", n_jobs=self.n_jobs, verbose=0) if g_tr is not None: grid.fit(X_tr[:, sorted_idx[:nfeat]], y_tr, groups=g_tr) else: grid.fit(X_tr[:, sorted_idx[:nfeat]], y_tr) curves[int(nfeat)] = float(grid.best_score_) if grid.best_score_ > best_cv: best_cv = float(grid.best_score_) best_n = int(nfeat) best_g = float(grid.best_params_["svm__gamma"]) return {"curves": curves, "best_nfeat": best_n, "best_cv": best_cv, "best_gamma": best_g} # ------------------------------------------------------------------ # Internals # ------------------------------------------------------------------ def _split(self, X, Y, groups, seed, split_type): """Image-level or patient-level 80:20 train/test split. Returns (X_tr, X_te, y_tr, y_te, g_tr, cv, n_tr_patients, n_te_patients); g_tr is None and cv is 5 for image-level splits. """ if split_type == "image": X_tr, X_te, y_tr, y_te = train_test_split( X, Y, test_size=0.2, random_state=seed, stratify=Y) return X_tr, X_te, y_tr, y_te, None, 5, None, None unique_g = np.unique(groups) g_labels = np.array([Y[groups == g][0] for g in unique_g]) tr_g, te_g = train_test_split( unique_g, test_size=0.2, random_state=seed, stratify=g_labels) tr_mask = np.isin(groups, tr_g) te_mask = np.isin(groups, te_g) cv = StratifiedGroupKFold(n_splits=5, shuffle=True, random_state=seed) return (X[tr_mask], X[te_mask], Y[tr_mask], Y[te_mask], groups[tr_mask], cv, len(np.unique(groups[tr_mask])), len(np.unique(groups[te_mask]))) def _rf_ranking(self, X_tr, y_tr): """Feature indices ranked by RF Gini importance (descending).""" fs = RandomForestClassifier(n_estimators=500, random_state=15, class_weight="balanced", n_jobs=self.n_jobs) fs.fit(X_tr, y_tr) return np.argsort(fs.feature_importances_)[::-1] def _grid_search(self, X_tr, y_tr, g_tr, cv, sorted_idx, nfeatures_list, gamma_logspace): """Search nfeatures × gamma; return (best_cv, best_nfeat, best_gamma).""" pipe = Pipeline([("scaler", StandardScaler()), ("svm", SVC(kernel="rbf", C=10))]) valid_nfeats = [n for n in nfeatures_list if n <= X_tr.shape[1]] best_cv, best_n, best_g = 0, None, None for nfeat in valid_nfeats: pg = {"svm__gamma": (1.0 / nfeat) * np.logspace(*gamma_logspace)} grid = GridSearchCV(pipe, pg, cv=cv, scoring="accuracy", n_jobs=self.n_jobs, verbose=0) if g_tr is not None: grid.fit(X_tr[:, sorted_idx[:nfeat]], y_tr, groups=g_tr) else: grid.fit(X_tr[:, sorted_idx[:nfeat]], y_tr) if grid.best_score_ > best_cv: best_cv = grid.best_score_ best_n = nfeat best_g = grid.best_params_["svm__gamma"] return best_cv, best_n, best_g def _load_data(self, model_name): data = np.load(os.path.join(self.features_dir, f"{model_name}_features.npz"), allow_pickle=True) X = data["X"] Y = data["Y"] filenames = data["filenames"] groups = np.array([self._image_to_patient.get(_f2n(f), f"unk_{i}") for i, f in enumerate(filenames)]) return X, Y, groups