1#!/usr/bin/env python32"""WP44 — semantic research mapping (RQ2/RQ5/RQ10 inputs): UMAP projection, HDBSCAN clusters, stability, temporal drift.3Inputs: data/processed/embeddings/<model>/shard_*.npy (+ ids), corpus_works, classification (primary categories).4Steps: (1) stratified sample of --n works (seed 44) for the 2-D atlas; (2) HDBSCAN on a 50-D UMAP of the sample, cluster labels5and top TF-IDF terms per cluster; (3) stability: ARI between HDBSCAN partitions of two bootstrap halves projected into the same6UMAP (repeated 5×); (4) category centroids per year in the native 768-D space, cosine drift year-over-year, inter-category7cosine similarity matrix per period (baseline vs post) → tests H4/H5 at the semantic level; (5) per-work novelty = 1 − max cosine8similarity to the previous 24 months' centroid of the same category (feeds experiments/emerging_detection.py).9Outputs: results/semantic/{atlas_sample.parquet, clusters.csv, cluster_terms.csv, stability.json, category_drift.csv,10category_similarity_<period>.csv, novelty.parquet}."""11import argparse, glob, json, os, numpy as np, pandas as pd, duckdb121314def load_embeddings(root, tag):15 E, I = [], []16 for f in sorted(glob.glob(f"{root}/data/processed/embeddings/{tag}/shard_*.npy")):17 if f.endswith(".ids.npy"):18 continue19 E.append(np.load(f)); I.append(np.load(f.replace(".npy", ".ids.npy")))20 return np.concatenate(E).astype(np.float32), np.concatenate(I)212223def main():24 ap = argparse.ArgumentParser(); ap.add_argument("--root", required=True); ap.add_argument("--model", default="allenai/specter2_base"); ap.add_argument("--n", type=int, default=200000)25 ap.add_argument("--seed", type=int, default=44); a = ap.parse_args(); R = a.root; tag = a.model.replace("/", "__"); out = f"{R}/results/semantic"; os.makedirs(out, exist_ok=True)26 E, ids = load_embeddings(R, tag); print("embeddings", E.shape)27 con = duckdb.connect(f"{R}/data/processed/wp44.duckdb", read_only=True) # read-only: other processes may hold the database28 meta = con.execute(f"SELECT c.id AS work_id, c.publication_year AS year, c.period, c.first_month, p.primary_category FROM corpus_works c LEFT JOIN read_parquet('{R}/data/processed/classification/primary/*.parquet') p ON p.work_id = c.id WHERE c.eligible").df().set_index("work_id")29 pos = {w: i for i, w in enumerate(ids)}; meta = meta[meta.index.isin(pos)]; meta["row"] = [pos[w] for w in meta.index]30 rng = np.random.default_rng(a.seed)31 # (4) category centroids per year + drift + similarity32 rows = []; cents = {}33 for (cat, y), g in meta.dropna(subset=["primary_category"]).groupby(["primary_category", "year"]):34 if len(g) < 100:35 continue36 v = E[g.row.values].mean(0); v /= np.linalg.norm(v) + 1e-9; cents[(cat, y)] = v37 for (cat, y), v in cents.items():38 prev = cents.get((cat, y - 1)); rows.append({"category": cat, "year": y, "n": int(((meta.primary_category == cat) & (meta.year == y)).sum()), "drift_cos_dist": None if prev is None else float(1 - v @ prev)})39 pd.DataFrame(rows).to_csv(f"{out}/category_drift.csv", index=False)40 for per, yrs in [("baseline", range(2018, 2023)), ("post", range(2023, 2027))]:41 cats = sorted({c for (c, y) in cents if y in yrs}); M = np.zeros((len(cats), len(cats)))42 for i, c1 in enumerate(cats):43 v1 = np.mean([cents[(c1, y)] for y in yrs if (c1, y) in cents], 0)44 for j, c2 in enumerate(cats):45 v2 = np.mean([cents[(c2, y)] for y in yrs if (c2, y) in cents], 0); M[i, j] = float(v1 @ v2 / (np.linalg.norm(v1) * np.linalg.norm(v2) + 1e-9))46 pd.DataFrame(M, index=cats, columns=cats).to_csv(f"{out}/category_similarity_{per}.csv")47 # (5) novelty vs trailing 24-month centroid of the same category48 meta_m = meta.dropna(subset=["primary_category"]).copy(); months = sorted(meta_m.first_month.unique()); nov = np.full(len(meta_m), np.nan, dtype=np.float32)49 mcent = {}50 for (cat, m), g in meta_m.groupby(["primary_category", "first_month"]):51 mcent[(cat, m)] = (E[g.row.values].sum(0), len(g))52 midx = {m: i for i, m in enumerate(months)}; arr_pos = {w: i for i, w in enumerate(meta_m.index)}53 for (cat, m), g in meta_m.groupby(["primary_category", "first_month"]):54 i0 = midx[m]; prev = [mcent[(cat, mm)] for mm in months[max(0, i0 - 24):i0] if (cat, mm) in mcent]55 if not prev:56 continue57 v = np.sum([p[0] for p in prev], 0) / max(sum(p[1] for p in prev), 1); v /= np.linalg.norm(v) + 1e-958 sims = E[g.row.values] @ v59 for w, s in zip(g.index, sims):60 nov[arr_pos[w]] = 1 - s61 meta_m["novelty"] = nov; meta_m[["novelty", "primary_category", "first_month"]].reset_index().to_parquet(f"{out}/novelty.parquet")62 # (1–3) atlas on a sample63 import umap, hdbscan64 from sklearn.metrics import adjusted_rand_score65 samp = meta.sample(min(a.n, len(meta)), random_state=a.seed); X = E[samp.row.values]66 red = umap.UMAP(n_components=50, n_neighbors=30, min_dist=0.0, metric="cosine", random_state=a.seed).fit_transform(X)67 lab = hdbscan.HDBSCAN(min_cluster_size=200, min_samples=20).fit_predict(red)68 xy = umap.UMAP(n_components=2, n_neighbors=30, min_dist=0.1, metric="cosine", random_state=a.seed).fit_transform(X)69 samp = samp.assign(cluster=lab, x=xy[:, 0], y=xy[:, 1]); samp.reset_index().to_parquet(f"{out}/atlas_sample.parquet")70 titles = con.execute(f"SELECT id AS work_id, title FROM corpus_works WHERE id IN (SELECT work_id FROM read_parquet('{out}/atlas_sample.parquet'))").df().set_index("work_id").title71 from sklearn.feature_extraction.text import TfidfVectorizer72 docs = samp.index.map(lambda w: titles.get(w, "") or ""); vec = TfidfVectorizer(max_features=50000, stop_words="english", ngram_range=(1, 2), min_df=5); T = vec.fit_transform(docs); terms = np.array(vec.get_feature_names_out())73 crow = []74 for c in sorted(set(lab)):75 m = lab == c; v = np.asarray(T[m].mean(0)).ravel(); top = terms[np.argsort(-v)[:12]]76 pc = samp[m].primary_category.value_counts(normalize=True)77 crow.append({"cluster": int(c), "n": int(m.sum()), "top_terms": "; ".join(top), "top_category": pc.index[0] if len(pc) else None, "top_category_share": float(pc.iloc[0]) if len(pc) else None,78 "share_post": float((samp[m].period == "post").mean())})79 pd.DataFrame(crow).to_csv(f"{out}/clusters.csv", index=False)80 aris = []81 for r in range(5):82 idx = rng.permutation(len(X)); h1, h2 = idx[: len(X) // 2], idx[len(X) // 2:]83 l1 = hdbscan.HDBSCAN(min_cluster_size=100, min_samples=20).fit_predict(red[h1]); l2 = hdbscan.HDBSCAN(min_cluster_size=100, min_samples=20).fit_predict(red[h2])84 # compare through nearest-neighbour label transfer from half 1 to half 285 from sklearn.neighbors import KNeighborsClassifier86 knn = KNeighborsClassifier(n_neighbors=5).fit(red[h1], l1); aris.append(float(adjusted_rand_score(l2, knn.predict(red[h2]))))87 json.dump({"n_sample": int(len(X)), "n_clusters": int(len(set(lab)) - (1 if -1 in lab else 0)), "noise_share": float((lab == -1).mean()), "ari_bootstrap": aris, "ari_mean": float(np.mean(aris))}, open(f"{out}/stability.json", "w"), indent=1)88 print("atlas done", len(set(lab)), "clusters; ARI", np.mean(aris))899091if __name__ == "__main__":92 main()93