{"cells":[{"cell_type":"markdown","id":"4a6caa1a","metadata":{},"source":"# RSNA 2026 Knee — a look at the data\n\nPoking around the RSNA 2026 knee data before committing to a pipeline.\n\n- Tables / metadata / report text: full pass (CSVs are cheap)\n- DICOM pixels: only a few series (pixel-level plots are marked *(sample)*)\n\nThe setup: 12-class multi-label knee MRI classification, metric is ROC AUC. The catch — out of 4407 studies only **58 have expert labels**; the rest come with just a multilingual radiology report. So the three things I actually wanted to understand:\n\n1. What does the little gold set look like (prevalence, co-occurrence, label structure)?\n2. How much signal can we squeeze out of 4407 free-text reports (language, length, duplication, keywords)?\n3. What are we feeding the network (series mix, geometry, contrast, scanner/protocol spread)?"},{"cell_type":"code","execution_count":1,"id":"081a7cc7","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:14.799507Z","iopub.status.busy":"2026-08-07T05:15:14.799028Z","iopub.status.idle":"2026-08-07T05:15:15.356909Z","shell.execute_reply":"2026-08-07T05:15:15.356698Z"}},"outputs":[],"source":"from pathlib import Path\nimport sys, re, math, unicodedata\nfrom collections import Counter\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nfrom matplotlib.gridspec import GridSpec\nfrom scipy import stats\n\n# --- path setup ---\n# On Kaggle the competition mounts at /kaggle/input/competitions/rsna-knee-abnormality-detection/.\n# Locally I keep the same layout (train.csv / *_series.csv / test_series/<study>/<series>/*.dcm),\n# so the notebook runs unchanged in both places.\nCANDIDATES = [\n    Path('/kaggle/input/competitions/rsna-knee-abnormality-detection'),\n    Path('/kaggle/input/rsna-knee-abnormality-detection'),\n    Path('../data'),\n    Path('/home/lian/notebooks/RSNA26/data'),\n]\nfor cand in CANDIDATES:\n    if (cand / 'train.csv').exists():\n        DATA = cand\n        break\nelse:\n    raise FileNotFoundError(f'train.csv not found under any of {CANDIDATES}')\nKAGGLE = str(DATA).startswith('/kaggle/')\nprint('DATA =', DATA if KAGGLE else DATA.resolve(), '| on Kaggle:', KAGGLE)\n\nLABELS = ['ACL','MCL','Medial Meniscus','Lateral Meniscus','Medial OA','Lateral OA',\n          'PF OA','Effusion','Synovitis',\"Baker's\",'Contusion','Fracture']\nN_CLASSES = len(LABELS)\n\npd.set_option('display.width', 200)\npd.set_option('display.max_columns', 40)\nplt.rcParams.update({'figure.dpi': 110, 'axes.grid': True, 'grid.alpha': 0.25,\n                     'axes.titlesize': 11, 'axes.labelsize': 10})\n\ndef cramers_v(confusion):\n    chi2 = stats.chi2_contingency(confusion, correction=False)[0]\n    n = confusion.sum()\n    return math.sqrt(max(0.0, chi2 / n))"},{"cell_type":"markdown","id":"f13af4b6","metadata":{},"source":"## 1. Study table: the gold set under a microscope"},{"cell_type":"code","execution_count":2,"id":"3386dd3b","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:15.358624Z","iopub.status.busy":"2026-08-07T05:15:15.358352Z","iopub.status.idle":"2026-08-07T05:15:15.402649Z","shell.execute_reply":"2026-08-07T05:15:15.402409Z"}},"outputs":[],"source":"train = pd.read_csv(DATA / 'train.csv')\nlabeled = train[train['ACL'].notna()].copy()\nY = labeled[LABELS].astype(int).values\n\nprint(f'studies          : {len(train)}')\nprint(f'expert-labeled   : {len(labeled)}  ({len(labeled)/len(train):.1%})')\nprint(f'report-only      : {len(train)-len(labeled)}')\nprint(f'report missing   : {train[\"Report\"].isna().sum()}')\nprint(f'duplicate UIDs   : {train[\"StudyInstanceUID\"].duplicated().sum()}')"},{"cell_type":"markdown","id":"5d6ea083","metadata":{},"source":"### 1.1 Prevalence with Wilson 95% CIs\n\nn=58, so point estimates wobble a lot. Wilson intervals behave better than the normal approx at counts like 1/58."},{"cell_type":"code","execution_count":3,"id":"56619d13","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:15.404003Z","iopub.status.busy":"2026-08-07T05:15:15.403915Z","iopub.status.idle":"2026-08-07T05:15:15.506135Z","shell.execute_reply":"2026-08-07T05:15:15.505767Z"}},"outputs":[],"source":"def wilson(k, n, z=1.959964):\n    p = k / n\n    denom = 1 + z**2 / n\n    center = (p + z**2 / (2*n)) / denom\n    half = z * math.sqrt(p*(1-p)/n + z**2/(4*n**2)) / denom\n    return center - half, center + half\n\nprev = pd.DataFrame({'k': labeled[LABELS].sum(), 'n': len(labeled)})\nprev['p'] = prev.k / prev.n\nci = prev.apply(lambda r: wilson(r.k, r.n), axis=1)\nprev['lo'] = [c[0] for c in ci]; prev['hi'] = [c[1] for c in ci]\nprev = prev.sort_values('p')\n\nfig, ax = plt.subplots(figsize=(9, 4.5))\ny = np.arange(len(prev))\nax.barh(y, prev.p * 100, alpha=0.75)\nax.errorbar(prev.p * 100, y,\n            xerr=[(prev.p - prev.lo) * 100, (prev.hi - prev.p) * 100],\n            fmt='none', ecolor='black', elinewidth=1.2, capsize=3)\nfor i, (p, k) in enumerate(zip(prev.p, prev.k)):\n    ax.text(p * 100 + 1, i, f'{p:.0%} (n={int(k)})', va='center', fontsize=8)\nax.set_yticks(y, prev.index)\nax.set_xlabel('prevalence %  (whiskers = Wilson 95% CI)')\nax.set_title(f'Label prevalence — {len(labeled)} gold studies')\nax.set_xlim(0, 105)\nplt.tight_layout(); plt.show()\n\nprint('Effective sample size vs. rarest class: MCL has', int(prev.loc['MCL','k']), 'positives.',\n      'A per-class AUC CI will be huge — aggregate metric hides this.')"},{"cell_type":"markdown","id":"67d301f5","metadata":{},"source":"### 1.2 Multi-label structure: cardinality, density, label-set entropy"},{"cell_type":"code","execution_count":4,"id":"a24e5d77","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:15.507749Z","iopub.status.busy":"2026-08-07T05:15:15.507631Z","iopub.status.idle":"2026-08-07T05:15:15.570149Z","shell.execute_reply":"2026-08-07T05:15:15.569954Z"}},"outputs":[],"source":"n_pos = Y.sum(axis=1)\ncard = n_pos.mean()\ndens = card / N_CLASSES\n\ncombos = Counter(map(tuple, Y.tolist()))\nent = stats.entropy(np.array(list(combos.values())))\n\nprint(f'label cardinality (mean positives/study): {card:.2f}')\nprint(f'label density     (cardinality / 12)    : {dens:.2f}')\nprint(f'distinct label sets: {len(combos)} / {len(labeled)} studies | entropy = {ent:.2f} nats')\nprint('most common combos:')\nfor combo, cnt in combos.most_common(5):\n    names = [LABELS[i] for i, b in enumerate(combo) if b]\n    print(f'  {cnt:3d}×  {names}')\n\nfig, ax = plt.subplots(figsize=(6.5, 4))\nax.hist(n_pos, bins=np.arange(-0.5, n_pos.max() + 1.5), rwidth=0.8)\nax.axvline(card, color='r', ls='--', lw=1, label=f'mean = {card:.2f}')\nax.set_xlabel('# positive labels'); ax.set_ylabel('# studies')\nax.set_title('Label cardinality distribution')\nax.legend()\nplt.tight_layout(); plt.show()"},{"cell_type":"markdown","id":"40c36905","metadata":{},"source":"### 1.3 Co-occurrence: counts, lift, Cramér's V\n\nRaw counts flatter pairs of common labels. **Lift** = P(A∩B)/P(A)P(B) (>1 = co-occur more than chance); **Cramér's V** gives a symmetric strength. Worth knowing which pairs the model gets almost for free via correlation."},{"cell_type":"code","execution_count":5,"id":"297ad573","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:15.571539Z","iopub.status.busy":"2026-08-07T05:15:15.571412Z","iopub.status.idle":"2026-08-07T05:15:15.849033Z","shell.execute_reply":"2026-08-07T05:15:15.84882Z"}},"outputs":[],"source":"n = len(labeled)\nco = labeled[LABELS].T.dot(labeled[LABELS]).values.astype(float)\np = Y.mean(axis=0)\nexp = np.outer(p, p) * n\nlift = co / np.maximum(exp, 1e-9)\n\nV = np.zeros((N_CLASSES, N_CLASSES))\nfor i in range(N_CLASSES):\n    for j in range(N_CLASSES):\n        a = Y[:, i]; b = Y[:, j]\n        tab = np.array([[((a==1)&(b==1)).sum(), ((a==1)&(b==0)).sum()],\n                        [((a==0)&(b==1)).sum(), ((a==0)&(b==0)).sum()]])\n        V[i, j] = cramers_v(tab)\n\nfig, axes = plt.subplots(1, 3, figsize=(17, 4.8))\nfor ax, M, title, cmap in [\n        (axes[0], co, 'co-occurrence count', 'viridis'),\n        (axes[1], lift, 'lift (diag = 1/p)', 'magma'),\n        (axes[2], V, \"Cramér's V\", 'cividis')]:\n    im = ax.imshow(M, cmap=cmap)\n    ax.set_xticks(range(N_CLASSES), LABELS, rotation=90, fontsize=7)\n    ax.set_yticks(range(N_CLASSES), LABELS, fontsize=7)\n    ax.set_title(title)\n    fig.colorbar(im, ax=ax, shrink=0.8)\nplt.tight_layout(); plt.show()\n\niu = np.triu_indices(N_CLASSES, 1)\ntop = sorted(zip(V[iu], lift[iu], co[iu], [LABELS[i] for i in iu[0]], [LABELS[j] for j in iu[1]]),\n             reverse=True)\nprint('Top associated pairs (V, lift, co-count):')\nfor v, l, c, a, b in top[:8]:\n    print(f'  V={v:.2f} lift={l:4.1f} n={int(c):2d}   {a} × {b}')"},{"cell_type":"markdown","id":"8a2836d3","metadata":{},"source":"### 1.4 Do we even have negatives? — per-class balance\n\nAUC needs both classes present. Checking nothing is all-positive/all-negative in the gold set, and how badly a grouped CV split starves rare classes."},{"cell_type":"code","execution_count":6,"id":"750d2d17","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:15.850681Z","iopub.status.busy":"2026-08-07T05:15:15.850529Z","iopub.status.idle":"2026-08-07T05:15:15.855503Z","shell.execute_reply":"2026-08-07T05:15:15.855267Z"}},"outputs":[],"source":"bal = pd.DataFrame({'pos': Y.sum(0), 'neg': (1 - Y).sum(0)}, index=LABELS)\nbal['pos_rate'] = bal.pos / len(labeled)\nbal['min_class'] = bal[['pos', 'neg']].min(axis=1)\nprint(bal.sort_values('min_class').to_string())\nprint()\nprint('Any class with <3 positives or <3 negatives:', (bal.min_class < 3).any())\nprint('In a 5-fold grouped CV (~11-12 studies/val), expected val positives for rarest class:',\n      f\"{bal.pos.min() / 5:.1f} — per-fold AUC for MCL will be noisy/meaningless; aggregate over folds or use all 58 for final.\")"},{"cell_type":"markdown","id":"8175e27e","metadata":{},"source":"## 2. Reports — the weak-label goldmine (4407 free-text, many languages)\n\nThe only supervision for 98.7% of studies. Profiling: language mix, length, duplication/template-ness, and crude per-class keyword hit rates."},{"cell_type":"code","execution_count":7,"id":"3f522505","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:15.856788Z","iopub.status.busy":"2026-08-07T05:15:15.856653Z","iopub.status.idle":"2026-08-07T05:15:16.087426Z","shell.execute_reply":"2026-08-07T05:15:16.087115Z"}},"outputs":[],"source":"rep = train['Report'].fillna('')\nrep_len = rep.str.len()\nrep_words = rep.str.split().str.len()\n\nprint('chars : mean %5.0f | p5 %4.0f | p50 %4.0f | p95 %5.0f | max %5d' % (\n    rep_len.mean(), rep_len.quantile(.05), rep_len.median(), rep_len.quantile(.95), rep_len.max()))\nprint('words : mean %5.0f | p50 %4.0f | max %4d' % (\n    rep_words.mean(), rep_words.median(), rep_words.max()))\n\nfig, axes = plt.subplots(1, 2, figsize=(13, 3.8))\naxes[0].hist(rep_len, bins=100)\naxes[0].set_xlabel('chars'); axes[0].set_title('Report length (chars)')\naxes[1].hist(np.log10(rep_len + 1), bins=60)\naxes[1].set_xlabel('log10(chars+1)'); axes[1].set_title('log-scale — exposes short-report modes')\nplt.tight_layout(); plt.show()"},{"cell_type":"markdown","id":"ec7b3bbb","metadata":{},"source":"### 2.1 Language distribution — heuristic guess\n\nCheap heuristic (script detection + stopword votes) — take it as indicative only. `de`/`la` collide across Romance languages; the real pipeline should use langdetect or an LLM."},{"cell_type":"code","execution_count":8,"id":"09c730b9","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:16.088991Z","iopub.status.busy":"2026-08-07T05:15:16.088849Z","iopub.status.idle":"2026-08-07T05:15:16.978981Z","shell.execute_reply":"2026-08-07T05:15:16.978553Z"}},"outputs":[],"source":"def guess_lang(s: str) -> str:\n    if not s:\n        return 'none'\n    probe = s[:3000]\n    if any('CYRILLIC' in unicodedata.name(c, '') for c in probe):\n        return 'ru/cyrillic'\n    if any('\\u4e00' <= c <= '\\u9fff' for c in probe):\n        return 'zh'\n    if any('\\u0600' <= c <= '\\u06ff' for c in probe):\n        return 'ar'\n    sl = ' ' + s.lower() + ' '\n    hits = {\n        'es': len(re.findall(r'\\b(rodilla|que|del|las|con|una|hay)\\b', sl)),\n        'de': len(re.findall(r'\\b(der|die|und|mit|des|knie|ein|ist)\\b', sl)),\n        'fr': len(re.findall(r'\\b(genou|les|des|une|avec|est|pas)\\b', sl)),\n        'tr': len(re.findall(r'\\b(diz|ve|ile|bir|için|olan|eklem)\\b', sl)),\n        'pt': len(re.findall(r'\\b(joelho|das|dos|com|uma|para|não)\\b', sl)),\n        'nl': len(re.findall(r'\\b(het|een|van|met|voor|knie|zijn)\\b', sl)),\n        'it': len(re.findall(r'\\b(ginocchio|della|che|con|per|una)\\b', sl)),\n        'en': len(re.findall(r'\\b(the|and|of|with|there|knee|joint)\\b', sl)),\n    }\n    best = max(hits, key=hits.get)\n    return best if hits[best] >= 3 else 'other'\n\nlang = rep.map(guess_lang)\nlang_counts = lang.value_counts()\nprint(lang_counts.to_string())\n\nfig, ax = plt.subplots(figsize=(7, 4))\nlang_counts.plot.barh(ax=ax)\nax.set_title('Report language (rough heuristic)')\nplt.tight_layout(); plt.show()"},{"cell_type":"markdown","id":"05b43346","metadata":{},"source":"### 2.2 Template-ness / duplication\n\nIf lots of reports are near-duplicates (same clinic template), the effective diversity of the weak-label corpus is well under 4407. Cheap proxy: exact dup rate + dup rate of the normalized first 200 chars."},{"cell_type":"code","execution_count":9,"id":"27fda9be","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:16.980555Z","iopub.status.busy":"2026-08-07T05:15:16.980443Z","iopub.status.idle":"2026-08-07T05:15:17.12094Z","shell.execute_reply":"2026-08-07T05:15:17.120651Z"}},"outputs":[],"source":"def normalize_txt(s):\n    s = unicodedata.normalize('NFKD', s).lower()\n    s = re.sub(r'\\[date\\]|\\[.*?\\]', ' ', s)\n    s = re.sub(r'\\W+', ' ', s)\n    return s.strip()\n\nexact_dup = rep.duplicated().sum()\nhead200 = rep.map(lambda s: normalize_txt(s)[:200])\nhead_dup = head200.duplicated().sum()\nuniq_heads = head200.nunique()\n\nprint(f'exact duplicate reports         : {exact_dup} ({exact_dup/len(rep):.1%})')\nprint(f'duplicate normalized 200-ch head: {head_dup} ({head_dup/len(rep):.1%})  | unique heads: {uniq_heads}')\n\ntop_heads = head200.value_counts().head(5)\nprint('\\nmost common report openings (normalized):')\nfor h, c in top_heads.items():\n    print(f'  {c:3d}×  {h[:90]!r}')"},{"cell_type":"markdown","id":"31cad0eb","metadata":{},"source":"### 2.3 Keyword probe — how much weak signal is even there?\n\nPer-class multilingual keyword lists. Hit rate = fraction of reports with at least one keyword for the class. **Not the label** — negation/context ignored — but it bounds what a keyword miner can see, and flags classes that barely show up in text."},{"cell_type":"code","execution_count":10,"id":"2c9c2031","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.122502Z","iopub.status.busy":"2026-08-07T05:15:17.122404Z","iopub.status.idle":"2026-08-07T05:15:17.320885Z","shell.execute_reply":"2026-08-07T05:15:17.320676Z"}},"outputs":[],"source":"KEYWORDS = {\n    'ACL':             ['acl', 'ligamento cruzado anterior', 'vorderes kreuzband', 'ligament croisé antérieur',\n                        'ön çapraz bağ', 'anterior cruciate'],\n    'MCL':             ['mcl', 'ligamento colateral medial', 'mediales kollateralband', 'ligament collatéral médial',\n                        'iç yan bağ', 'medial collateral'],\n    'Medial Meniscus': ['menisco medial', 'menisco interno', 'innenmeniskus', 'ménisque médial',\n                        'medial menisküs', 'medial meniscus'],\n    'Lateral Meniscus':['menisco lateral', 'außenmeniskus', 'ménisque latéral', 'lateral menisküs',\n                        'lateral meniscus'],\n    'Medial OA':       ['artrosis femorotibial medial', 'gonarthrose medial', 'medial osteoarthritis',\n                        'medial compartment osteoarthritis'],\n    'Lateral OA':      ['artrosis femorotibial lateral', 'gonarthrose lateral', 'lateral osteoarthritis'],\n    'PF OA':           ['femoropatelar', 'patellofemoral', 'retropatellararthrose', 'fémoropatellaire'],\n    'Effusion':        ['derrame', 'effusion', 'erguss', 'épanchement', 'efüzyon', 'joint fluid'],\n    'Synovitis':       ['sinovitis', 'synovitis', 'synovialitis', 'synovite', 'sinovit'],\n    \"Baker's\":        ['baker', 'popliteal cyst', 'quiste de baker', 'bakersche zyste', 'kyste de baker'],\n    'Contusion':       ['contusion', 'contusión', 'prellung', 'bone bruise', 'ödeme óseo', 'edema oseo'],\n    'Fracture':        ['fractura', 'fraktur', 'fracture', 'kırık', 'broken'],\n}\n\nrep_l = rep.str.lower()\nhits = pd.DataFrame({c: rep_l.map(lambda s, kws=kws: any(k in s for k in kws))\n                     for c, kws in KEYWORDS.items()})\n\nhit_rate = hits.mean().sort_values()\nfig, ax = plt.subplots(figsize=(9, 4.5))\nax.barh(hit_rate.index, hit_rate.values * 100)\nfor i, v in enumerate(hit_rate.values * 100):\n    ax.text(v + 0.3, i, f'{v:.0f}%', va='center', fontsize=8)\nax.set_xlabel('% of reports with ≥1 class keyword (context-blind)')\nax.set_title('Keyword visibility per class — upper bound for a keyword weak-labeler')\nplt.tight_layout(); plt.show()\n\nprint('classes with <5% keyword visibility — these need synonyms/LLM mining or will be invisible:')\nprint(list(hit_rate[hit_rate < 0.05].index))"},{"cell_type":"markdown","id":"2ebb29ba","metadata":{},"source":"### 2.4 Keyword-vs-gold sanity check\n\nOn the 58 gold studies we can actually score the keyword probe against truth. n is tiny so it's noisy, but a keyword set whose recall is ~0 on gold is definitely broken."},{"cell_type":"code","execution_count":11,"id":"f8968371","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.322454Z","iopub.status.busy":"2026-08-07T05:15:17.32234Z","iopub.status.idle":"2026-08-07T05:15:17.32731Z","shell.execute_reply":"2026-08-07T05:15:17.327101Z"}},"outputs":[],"source":"gold_hits = hits.loc[labeled.index].values\nrows = []\nfor i, c in enumerate(LABELS):\n    y = Y[:, i]; h = gold_hits[:, i]\n    tp = int(((y==1)&(h==1)).sum()); fp = int(((y==0)&(h==1)).sum())\n    fn = int(((y==1)&(h==0)).sum()); tn = int(((y==0)&(h==0)).sum())\n    rec = tp / max(tp+fn, 1); prec = tp / max(tp+fp, 1)\n    rows.append({'class': c, 'gold_pos': int(y.sum()), 'kw_hit': int(h.sum()),\n                 'kw_recall': round(rec, 2), 'kw_prec': round(prec, 2),\n                 'tp': tp, 'fp': fp, 'fn': fn})\nkw_eval = pd.DataFrame(rows).set_index('class')\nprint(kw_eval.to_string())\nprint()\nprint('Caveat: n=58, and keyword lists are v1 — low recall mostly means missing synonyms,',\n      'and low precision often reflects negation (\"no effusion\") we ignore here.')"},{"cell_type":"markdown","id":"afd6b393","metadata":{},"source":"## 3. Series metadata — what the model actually sees\n\n24,371 series across 4407 studies. Each study is a bag of series; the model has to aggregate across planes/weightings."},{"cell_type":"code","execution_count":12,"id":"2df4d3d8","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.328575Z","iopub.status.busy":"2026-08-07T05:15:17.328392Z","iopub.status.idle":"2026-08-07T05:15:17.517998Z","shell.execute_reply":"2026-08-07T05:15:17.517767Z"}},"outputs":[],"source":"ts = pd.read_csv(DATA / 'train_series.csv')\nper_study = ts.groupby('StudyInstanceUID').size()\n\nprint(f'series rows          : {len(ts)}')\nprint(f'series/study         : mean {per_study.mean():.2f} | min {per_study.min()} | max {per_study.max()}')\nprint(f'Fluid==FatSuppression: {(ts.Fluid_Sensitive == ts.Fat_Suppression).all()}  -> single \"sequence-type\" flag')\n\nfig, axes = plt.subplots(1, 3, figsize=(16, 4))\n\naxes[0].hist(per_study, bins=np.arange(0.5, per_study.max()+1.5), rwidth=0.8)\naxes[0].set_xlabel('# series in study'); axes[0].set_title('Series per study')\n\npd.crosstab(ts.Anatomical_Plane, ts.Fluid_Sensitive).plot.bar(ax=axes[1])\naxes[1].set_title('Plane × Fluid_Sensitive'); axes[1].set_xlabel('')\n\nts['stype'] = ts.Anatomical_Plane.str[:3] + np.where(ts.Fluid_Sensitive==1, '+FS', '-FS')\ncover = ts.groupby('StudyInstanceUID')['stype'].nunique()\naxes[2].hist(cover, bins=np.arange(0.5, cover.max()+1.5), rwidth=0.8)\naxes[2].set_xlabel('# distinct series-types'); axes[2].set_title('Series-type diversity per study')\nplt.tight_layout(); plt.show()\n\nprint('\\nSeries-type mix (plane × fluid-sensitive):')\nprint(ts.stype.value_counts().to_string())"},{"cell_type":"code","execution_count":13,"id":"a5bbf2d0","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.519494Z","iopub.status.busy":"2026-08-07T05:15:17.519407Z","iopub.status.idle":"2026-08-07T05:15:17.553861Z","shell.execute_reply":"2026-08-07T05:15:17.553667Z"}},"outputs":[],"source":"pivot = ts.pivot_table(index='Anatomical_Plane', columns='Fluid_Sensitive',\n                       values='SeriesInstanceUID', aggfunc='count').fillna(0).astype(int)\nprint('Series count by plane × fluid:')\nprint(pivot)\n\nplanes_per_study = ts.groupby('StudyInstanceUID')['Anatomical_Plane'].agg(lambda s: set(s))\nhave_all3 = planes_per_study.map(lambda s: {'Sagittal','Coronal','Axial'} <= set(s)).mean()\nprint(f'\\nstudies with all 3 planes: {have_all3:.1%}')\nprint('plane-set distribution:')\nprint(planes_per_study.map(lambda s: '+'.join(sorted(x[:3] for x in s))).value_counts().head(8).to_string())"},{"cell_type":"markdown","id":"781bf50a","metadata":{},"source":"## 4. DICOM pixels & headers — a few series, in detail\n\nDICOM helpers inlined below so this runs standalone on Kaggle. We look at:\n\n- **Geometry sanity**: does `Anatomical_Plane` agree with `ImageOrientationPatient`? (spoiler: yes, and the slices are slightly oblique)\n- **Acquisition params**: TR/TE/ETL/flip — can we read weighting off the headers?\n- **Foreground/background**: how much of the FOV is air (zeros)? — relevant to the efficiency prize\n- **Slice spacing**: thickness vs. spacing (the gap)"},{"cell_type":"code","execution_count":14,"id":"a0af8991","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.555413Z","iopub.status.busy":"2026-08-07T05:15:17.555301Z","iopub.status.idle":"2026-08-07T05:15:17.652389Z","shell.execute_reply":"2026-08-07T05:15:17.651963Z"}},"outputs":[],"source":"import pydicom\nfrom pydicom.pixel_data_handlers.util import apply_voi_lut\n\ndef load_dicom_pixel(path):\n    \"\"\"pixel array (float32) + dataset, with VOI LUT + rescale applied.\"\"\"\n    ds = pydicom.dcmread(str(path))\n    arr = ds.pixel_array\n    try:\n        arr = apply_voi_lut(arr, ds)\n    except Exception:\n        pass\n    arr = arr.astype(np.float32)\n    slope = float(getattr(ds, 'RescaleSlope', 1.0) or 1.0)\n    intercept = float(getattr(ds, 'RescaleIntercept', 0.0) or 0.0)\n    return arr * slope + intercept, ds\n\ndef series_files(series_dir):\n    p = Path(series_dir)\n    return sorted([f for f in p.iterdir() if f.suffix.lower() == '.dcm'])\n\ndef load_series_volume(series_dir):\n    \"\"\"stack slices into (D,H,W) float32, sorted by position on slice normal.\"\"\"\n    files = series_files(series_dir)\n    if not files:\n        return None, []\n    slices, keys, dss = [], [], []\n    for f in files:\n        try:\n            arr, ds = load_dicom_pixel(f)\n        except Exception:\n            continue\n        if arr.ndim == 3:\n            arr = arr[0] if arr.shape[0] < arr.shape[-1] else arr[..., 0]\n        slices.append(arr)\n        ipp = getattr(ds, 'ImagePositionPatient', None)\n        iop = getattr(ds, 'ImageOrientationPatient', None)\n        if ipp is not None and iop is not None:\n            iop = np.asarray([float(x) for x in iop])\n            normal = np.cross(iop[:3], iop[3:])\n            keys.append(float(np.dot(np.asarray([float(x) for x in ipp]), normal)))\n        else:\n            keys.append(float(getattr(ds, 'InstanceNumber', 0) or 0))\n        dss.append(ds)\n    if not slices:\n        return None, []\n    order = np.argsort(keys)\n    vol = np.stack([slices[i] for i in order], axis=0).astype(np.float32)\n    return vol, [dss[i] for i in order]\n\ndef robust_normalize(vol, lo=0.5, hi=99.5):\n    \"\"\"percentile-clip + scale to [0,1] — robust across vendors/field strengths.\"\"\"\n    v = vol[np.isfinite(vol)]\n    if v.size == 0:\n        return np.zeros_like(vol, dtype=np.float32)\n    a, b = np.percentile(v, [lo, hi])\n    if b <= a:\n        b = a + 1e-6\n    return np.clip((vol - a) / (b - a), 0, 1).astype(np.float32)\n\nprint('dicom helpers ready')"},{"cell_type":"code","execution_count":15,"id":"dae12925","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.6538Z","iopub.status.busy":"2026-08-07T05:15:17.653716Z","iopub.status.idle":"2026-08-07T05:15:17.657534Z","shell.execute_reply":"2026-08-07T05:15:17.657327Z"}},"outputs":[],"source":"# find a few series of DICOMs. The layout is <DATA>/test_series/<study>/<series>/*.dcm both\n# on Kaggle and locally, so the same walk works — on the full train set we'd point this at train_series/.\ndef find_series_dirs(root, min_files=5, limit=3):\n    hits = {}\n    for f in Path(root).rglob('*.dcm'):\n        hits[f.parent] = hits.get(f.parent, 0) + 1\n    dirs = [d for d, c in hits.items() if c >= min_files]\n    return sorted(dirs)[:limit]\n\nseries_dirs = find_series_dirs(DATA / 'test_series', limit=3)\n# fall back to train_series/ if test_series is empty (e.g. full data locally)\nif not series_dirs and (DATA / 'train_series').exists():\n    series_dirs = find_series_dirs(DATA / 'train_series', limit=3)\n\nprint(f'{len(series_dirs)} series to inspect')\nfor d in series_dirs:\n    print(' ', d.relative_to(DATA), f'({len(list(d.glob(\"*.dcm\")))} dicoms)')"},{"cell_type":"code","execution_count":16,"id":"ee610305","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:17.658711Z","iopub.status.busy":"2026-08-07T05:15:17.658625Z","iopub.status.idle":"2026-08-07T05:15:18.786877Z","shell.execute_reply":"2026-08-07T05:15:18.786454Z"}},"outputs":[],"source":"def plane_from_iop(iop):\n    \"\"\"geometric plane from row/col direction cosines: which axis the slice normal follows.\"\"\"\n    iop = np.asarray([float(x) for x in iop])\n    normal = np.cross(iop[:3], iop[3:])\n    idx = int(np.argmax(np.abs(normal)))\n    return ['Sagittal', 'Coronal', 'Axial'][idx], normal\n\nrows = []\nvolumes = {}\nfor d in series_dirs:\n    vol, dss = load_series_volume(d)\n    if vol is None: continue\n    name = str(d)\n    volumes[name] = vol\n    ds0 = dss[0]\n    iop = getattr(ds0, 'ImageOrientationPatient', None)\n    geom, normal = plane_from_iop(iop) if iop is not None else ('?', None)\n    thick = float(getattr(ds0, 'SliceThickness', np.nan) or np.nan)\n    space = float(getattr(ds0, 'SpacingBetweenSlices', np.nan) or np.nan)\n    rows.append({\n        'series': name[-50:], 'slices': vol.shape[0], 'H': vol.shape[1], 'W': vol.shape[2],\n        'SeriesDesc': str(getattr(ds0, 'SeriesDescription', '')),\n        'geom_plane': geom,\n        'normal': None if normal is None else np.round(normal, 2).tolist(),\n        'pix_mm': [round(float(x),4) for x in getattr(ds0, 'PixelSpacing', [np.nan,np.nan])],\n        'thick_mm': thick, 'space_mm': space, 'gap_mm': round(space-thick, 2),\n        'TR': float(getattr(ds0, 'RepetitionTime', np.nan) or np.nan),\n        'TE': float(getattr(ds0, 'EchoTime', np.nan) or np.nan),\n        'ETL': int(getattr(ds0, 'EchoTrainLength', 0) or 0),\n        'flip': float(getattr(ds0, 'FlipAngle', np.nan) or np.nan),\n        'field_T': float(getattr(ds0, 'MagneticFieldStrength', np.nan) or np.nan),\n        'model': str(getattr(ds0, 'ManufacturerModelName', '')),\n        'coil': str(getattr(ds0, 'ReceiveCoilName', '')),\n        'bg_frac': round(float((vol == 0).mean()), 3),\n        'p99.9': int(np.percentile(vol, 99.9)), 'max': int(vol.max()),\n    })\n\nmeta = pd.DataFrame(rows)\nwith pd.option_context('display.max_colwidth', 55):\n    display(meta)"},{"cell_type":"markdown","id":"03ee759c","metadata":{},"source":"### 4.1 Slice mosaics + central-slice detail\n\nNote the big black FOV margins and (in sagittal views) how little of the frame is knee. An ROI crop is basically free accuracy-per-FLOP."},{"cell_type":"code","execution_count":17,"id":"82425003","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:18.788352Z","iopub.status.busy":"2026-08-07T05:15:18.788262Z","iopub.status.idle":"2026-08-07T05:15:20.294032Z","shell.execute_reply":"2026-08-07T05:15:20.293436Z"}},"outputs":[],"source":"for name, vol in volumes.items():\n    vn = robust_normalize(vol)\n    D = vn.shape[0]\n    idx = np.linspace(0, D-1, min(10, D)).round().astype(int)\n    fig, axes = plt.subplots(1, len(idx), figsize=(17, 2.1))\n    for ax, i in zip(np.atleast_1d(axes), idx):\n        ax.imshow(vn[i], cmap='gray', interpolation='nearest')\n        ax.set_title(f'#{i}', fontsize=8); ax.axis('off')\n    fig.suptitle(f'{name[-60:]}   ({vol.shape[0]}×{vol.shape[1]}×{vol.shape[2]})', fontsize=9)\n    plt.tight_layout(); plt.show()"},{"cell_type":"code","execution_count":18,"id":"65f1f977","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:20.295606Z","iopub.status.busy":"2026-08-07T05:15:20.295481Z","iopub.status.idle":"2026-08-07T05:15:21.764293Z","shell.execute_reply":"2026-08-07T05:15:21.764047Z"}},"outputs":[],"source":"fig = plt.figure(figsize=(15, 4.2 * len(volumes)))\ngs = GridSpec(len(volumes), 3, figure=fig)\n\nfor r, (name, vol) in enumerate(volumes.items()):\n    mid = vol.shape[0] // 2\n    raw = vol[mid]\n    norm = robust_normalize(vol)[mid]\n    mask = raw > 0.01 * raw.max()\n    ys, xs = np.where(mask)\n    bbox = (xs.min(), ys.min(), xs.max(), ys.max()) if len(xs) else (0,0,raw.shape[1]-1, raw.shape[0]-1)\n\n    ax0 = fig.add_subplot(gs[r, 0]); ax0.imshow(raw, cmap='gray'); ax0.set_title(f'raw #{mid}')\n    ax1 = fig.add_subplot(gs[r, 1]); ax1.imshow(norm, cmap='gray'); ax1.set_title('percentile-norm')\n    ax2 = fig.add_subplot(gs[r, 2])\n    x0,y0,x1,y1 = bbox\n    ax2.imshow(norm[y0:y1+1, x0:x1+1], cmap='gray'); ax2.set_title(f'FG crop {x1-x0+1}×{y1-y0+1}')\n    for ax in (ax0, ax1, ax2): ax.axis('off')\n    fg = mask.mean()\n    ax0.set_ylabel(name.split('/')[-1][:24], fontsize=7, rotation=0, labelpad=60, va='center')\n    print(f'{name[-50:]}  foreground={fg:.1%}  bbox={bbox}  crop-side≈{max(x1-x0, y1-y0)}px')\nplt.tight_layout(); plt.show()"},{"cell_type":"markdown","id":"a3cb2efb","metadata":{},"source":"### 4.2 Intensity distributions & normalization\n\nRaw uint16 ranges differ per series; percentile normalization equalizes scale, but watch the **big zero-background mass** — include background in the percentiles and you waste dynamic range on air."},{"cell_type":"code","execution_count":19,"id":"b3d5a861","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:21.768859Z","iopub.status.busy":"2026-08-07T05:15:21.768756Z","iopub.status.idle":"2026-08-07T05:15:22.196261Z","shell.execute_reply":"2026-08-07T05:15:22.195885Z"}},"outputs":[],"source":"fig, axes = plt.subplots(1, 2, figsize=(14, 4))\nfor name, vol in volumes.items():\n    v = vol.ravel()\n    step = max(1, v.size // 200000)\n    axes[0].hist(v[::step], bins=120, histtype='step', density=True, label=name[-38:])\n    fg = v[v > 0]\n    axes[1].hist(fg[::step], bins=120, histtype='step', density=True, label=name[-38:])\naxes[0].set_title('all voxels (note bg spike)'); axes[0].set_yscale('log')\naxes[1].set_title('foreground only'); axes[1].set_yscale('log')\nfor ax in axes:\n    ax.set_xlabel('raw pixel value')\n    h, l = ax.get_legend_handles_labels()\n    if h: ax.legend(fontsize=6)\nplt.tight_layout(); plt.show()"},{"cell_type":"markdown","id":"47cf43bc","metadata":{},"source":"### 4.3 Geometry check: does `Anatomical_Plane` match `ImageOrientationPatient`?\n\nInferring the geometric plane from the slice normal and comparing with the CSV label. (On Kaggle we match series by UID against `train_series.csv`/`test_series.csv`.)"},{"cell_type":"code","execution_count":20,"id":"06a0f510","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:22.197719Z","iopub.status.busy":"2026-08-07T05:15:22.197633Z","iopub.status.idle":"2026-08-07T05:15:22.868449Z","shell.execute_reply":"2026-08-07T05:15:22.868013Z"}},"outputs":[],"source":"test_series = pd.read_csv(DATA / 'test_series.csv')\n# on Kaggle the sample dicoms may belong to train or test; join both to look up plane labels\nall_series = pd.concat([\n    ts.assign(split='train'),\n    test_series.assign(split='test'),\n], ignore_index=True)\n\nchecks = []\nfor d in series_dirs:\n    suid = d.name\n    row = all_series[all_series.SeriesInstanceUID == suid]\n    vol, dss = load_series_volume(d)\n    geom, normal = plane_from_iop(getattr(dss[0], 'ImageOrientationPatient'))\n    checks.append({\n        'series_uid': suid[-24:],\n        'csv_plane': row.Anatomical_Plane.iloc[0] if len(row) else '?',\n        'geom_plane': geom,\n        'normal_xyz': np.round(normal, 2).tolist(),\n        'match': bool(len(row) and row.Anatomical_Plane.iloc[0] == geom),\n        'SeriesDesc': str(getattr(dss[0], 'SeriesDescription', '')),\n    })\nchk = pd.DataFrame(checks)\nprint(chk.to_string(index=False))\nprint()\nprint('obliqueness: normals are not axis-aligned — the scanner prescribed slightly off-axis.',\n      'Fine for 2.5D CNNs; matters if you ever resample to a canonical grid.')"},{"cell_type":"markdown","id":"ffd7fcb2","metadata":{},"source":"## 5. Test set & submission format"},{"cell_type":"code","execution_count":21,"id":"70f1e8f4","metadata":{"execution":{"iopub.execute_input":"2026-08-07T05:15:22.869955Z","iopub.status.busy":"2026-08-07T05:15:22.869862Z","iopub.status.idle":"2026-08-07T05:15:22.877758Z","shell.execute_reply":"2026-08-07T05:15:22.877573Z"}},"outputs":[],"source":"test = pd.read_csv(DATA / 'test.csv')\nsub = pd.read_csv(DATA / 'sample_submission.csv')\n\nprint('test studies:', len(test), '| test series:', len(test_series))\nprint('submission: one row per StudyInstanceUID, 12 probability columns')\nwith pd.option_context('display.max_columns', 15):\n    display(sub.head(3))\n\nprint('test series per study:')\nprint(test_series.groupby('StudyInstanceUID').size().to_string())"},{"cell_type":"markdown","id":"e891eb08","metadata":{},"source":"## 6. What I'd actually take away from this\n\n**The evaluation is the biggest trap.**\nn=58, and the rarest class (MCL) has 9 positives. A 5-fold CV puts ~2 MCL positives in a val fold, so per-fold per-class AUC is basically noise, and every prevalence has a Wilson CI ±10–15 points wide. Pool predictions across folds before computing AUC, report CIs, and don't read anything into a single-split 0.02 gain.\n\n**The keyword weak-labeler (v1) is broken for exactly the classes we care about.**\nScored against the 58 gold studies: all three OA subtypes (Medial/Lateral/PF) get recall = 0 — the multilingual word lists just don't include the phrases radiologists use (e.g. Spanish \"artrosis\"). Effusion/meniscus are text-visible (recall ~0.6–0.7) but precision is capped by negation we don't handle (\"no effusion\"). So: keyword mining is at best a positive-precision seed; the real lever is LLM distillation with negation/section awareness. Now there's numbers behind that instead of a hunch.\n\n**Label structure is exploitable, with a caveat.**\n~4.1 positives/study, density 0.34, and co-occurrence follows anatomy: meniscus↔same-side OA, Effusion×Synovitis×Baker's (inflammation triad), ACL×Contusion (trauma). A shared multi-label representation should help; the flip side is errors correlate too, so one broken class can drag the aggregate AUC.\n\n**The input side is tidy, and there's a free FLOP win.**\nEvery study has all 3 planes (100%), `Fluid_Sensitive == Fat_Suppression` always, ~5.5 series/study → one sequence-type flag and a straightforward bag-of-series model. Sample DICOMs: Siemens 1.5T, 15-ch knee coil, ~960×960 uint16 (12-bit stored), ~30 slices, 3.0/3.3mm spacing, TR/TE consistent with PD/T2 TSE. Slices are a few degrees oblique — irrelevant for a 2.5D CNN, only matters if resampling. **~15–20% of voxels are pure background and the knee fills well under half the FOV** → a foreground crop is the cheapest accuracy-per-FLOP available for the efficiency prize.\n\n**The report corpus is smaller than it looks.**\nMedian ~1k chars, but 3% exact dups and ~10% share a normalized 200-char opening (biggest template cluster appears 37×). Effective diversity < 4407, spread thin across es/de/tr/en/fr/ru/… Quality of mining (LLM, negation) beats chasing coverage.\n\nSo the to-do list more or less writes itself:\n1. Fix the weak-label miner — synonym expansion + LLM distillation + negation, OA subtypes first (that's the blind spot).\n2. Rebuild evaluation around pooled-AUC + bootstrap CIs; stop trusting single folds.\n3. Add a foreground bbox crop to the DICOM pipeline (percentile normalization is already there) — nearly free FLOPs.\n4. Run the header/pixel scan over the full DICOM set to confirm the plane labels and background fraction hold outside the sample."}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.12.2"}},"nbformat":4,"nbformat_minor":5}