{"cells":[{"cell_type":"markdown","metadata":{},"source":"# 膝蓋 MRI 異常辨識｜正式提交版本\n## 載入已訓練模型，對 Kaggle 隱藏測試資料進行推論\n\n這一份只用於提交比賽。我先在右側 **Add Input** 加入訓練 Notebook `colachiang/mri-test` 的最新完整 Output，並保留官方比賽資料。接著關閉 Internet，選擇 GPU，使用 **Save Version → Save & Run All**。\n\n程式會自動尋找 Output 中的 `knee_manifest.json`，載入 ResNet-18 與五個 attention fold 模型。正式提交時，Kaggle 會把三筆示例測試資料換成約 1,300 筆隱藏測試資料，最後在 `/kaggle/working/submission.csv` 產生提交檔。若程式顯示找到 0 或多於 1 份模型，請只保留一份訓練 Output，或在下一格填入 `MODEL_DIR`。\n","id":"knee-000"},{"cell_type":"markdown","metadata":{},"source":"## 01｜先確定我要解決的問題\n\n官方任務是以一次 MRI 檢查（study）為單位，輸出 12 種異常的分數。同一個人可能同時有不只一種異常，所以我使用多標籤分類，最後用 sigmoid，不能用互斥的 softmax。\n\n| 英文欄位 | 我理解的辨識目標 |\n|---|---|\n| ACL / MCL | 前十字韌帶／內側副韌帶損傷 |\n| Medial Meniscus / Lateral Meniscus | 內側／外側半月板撕裂 |\n| Medial OA / Lateral OA / PF OA | 內側／外側脛股關節、髕股關節退化 |\n| Effusion / Synovitis | 關節積液／滑膜炎 |\n| Baker's | 貝克氏囊腫 |\n| Contusion / Fracture | 骨挫傷／骨折 |\n\n官方只提供少量人工標籤，其餘檢查附有多語言報告。**缺標籤不代表正常。** 第一版先使用有人工標籤的樣本完成影像基準；報告在這版只用於資料探索，不輸入影像模型，也不自動改寫成人工標籤。因此這版是影像模型，不宣稱已完成多模態訓練。\n\n正式測試沒有報告。評分是 12 類 ROC-AUC 的平均，並非 accuracy。官方資料頁約 570 GB，我直接讀取 Kaggle 掛載資料，不下載整包、不列印所有 DICOM 路徑。\n","id":"knee-001"},{"cell_type":"markdown","metadata":{},"source":"## 02｜我的執行方式\n\n第一次：我先在比賽頁加入比賽並接受規則，建立 Notebook，匯入這份 `.ipynb`，確認右側 Input 有本比賽資料。選擇 **GPU T4 x2**，訓練時開啟 Internet，保持 `MODE=\"train\"`，從頭執行。這版主要使用第一張 GPU；有兩張卡不代表自動加倍速度。\n\n訓練完成後：我用 **Save Version → Save & Run All** 保存輸出。`/kaggle/working` 中的即時 checkpoint 只保護仍存在的工作環境，**不是跨工作階段的雲端備份**；我會確認已完成的版本有保留輸出。\n\n正式提交：我複製這份 Notebook，將訓練版本的 Output 加入 Input，設 `MODE=\"infer\"`，關閉 Internet，再 Save & Run All 並提交 `submission.csv`。推論版不會重新訓練，也不依賴 train.csv。若找到多份模型，我會明確填入 `MODEL_DIR`，不讓程式猜。\n\n中斷續跑：同一工作階段直接重跑；新工作階段則掛載上一個已保存的 Output，設 `RESUME_DIR` 為其 `knee_project` 目錄。相同參數下會重用特徵與每一 epoch 的狀態。修改實驗設定時，改用新的輸出目錄，不混用舊結果。\n\n比賽官方頁目前列出：報名期限 2026-10-15 23:59 UTC、提交期限 2026-10-22 23:59 UTC、Notebook 執行上限 9 小時。提交前我仍以官方最新公告為準。\n","id":"knee-002"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"from pathlib import Path\nimport os, sys, json, time, random, hashlib, shutil, subprocess, importlib.util\n\nMODE = \"infer\"                    # 第二次正式提交改成 \"infer\"\nMODEL_DIR = \"\"                    # 推論：附加輸出的 knee_project 絕對路徑；只有一份時自動找\nRESUME_DIR = \"\"                   # 續跑：上一份輸出的 knee_project 絕對路徑\nOUT = Path(\"/kaggle/working/knee_project\")\nCFG = dict(seed=42, image_size=224, slices_per_plane=8, folds=5,\n           epochs=40, batch_size=8, lr=0.0003, weight_decay=0.01,\n           dropout=0.35, feature_batch_size=16, bootstrap=300,\n           model_version=\"knee-mil-v1\")\nTARGETS = [\"ACL\", \"MCL\", \"Medial Meniscus\", \"Lateral Meniscus\", \"Medial OA\",\n           \"Lateral OA\", \"PF OA\", \"Effusion\", \"Synovitis\", \"Baker's\", \"Contusion\", \"Fracture\"]\nPLANES = [\"Sagittal\", \"Coronal\", \"Axial\"]\nIS_SCORING = os.environ.get(\"KAGGLE_IS_COMPETITION_RERUN\", \"\").lower() in {\"1\", \"true\"}\nif IS_SCORING:\n    MODE = \"infer\"\nassert MODE in {\"train\", \"infer\"}\nINPUT = Path(\"/kaggle/input\")\nOUT.mkdir(parents=True, exist_ok=True)\nSOURCE = None\nif MODE == \"infer\":\n    candidates = ([Path(MODEL_DIR)] if MODEL_DIR else\n                  sorted({p.parent for p in INPUT.rglob(\"knee_manifest.json\")}))\n    if len(candidates) != 1:\n        raise RuntimeError(f\"需要唯一一份訓練輸出，找到 {len(candidates)} 份。請掛載 Output 並填 MODEL_DIR。\")\n    SOURCE = candidates[0]\n    manifest = json.loads((SOURCE / \"knee_manifest.json\").read_text())\n    assert manifest[\"targets\"] == TARGETS\n    CFG = manifest[\"config\"]\nelif RESUME_DIR:\n    SOURCE = Path(RESUME_DIR)\n\nprint(\"執行模式：\", MODE, \"｜官方隱藏測試：\", IS_SCORING)\nprint(\"輸出位置：\", OUT)\nprint(json.dumps(CFG, indent=2, ensure_ascii=False))\n","id":"knee-003"},{"cell_type":"markdown","metadata":{},"source":"## 03｜環境與離線依賴\n\n影像含有 JPEG Lossless、JPEG 2000 等壓縮格式，我使用 pydicom 搭配 GDCM 解碼。訓練時會保存相同版本的 wheel；提交時只從掛載的 wheel 安裝，不連線。PyTorch、torchvision、NumPy 與 scikit-learn 使用 Kaggle 既有環境。若未來 Kaggle 更換 Python 版本，我需要重新產生對應的 wheel。\n","id":"knee-004"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"import importlib.metadata as metadata\n\npackages = {\"pydicom\": \"pydicom\", \"gdcm\": \"python-gdcm\"}\nmissing = [dist for module, dist in packages.items() if importlib.util.find_spec(module) is None]\nif missing:\n    args = [sys.executable, \"-m\", \"pip\", \"install\", \"--quiet\", \"--no-deps\"]\n    if MODE == \"infer\":\n        args += [\"--no-index\", \"--find-links\", str(SOURCE / \"wheels\")]\n    subprocess.check_call(args + missing)\n\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport pydicom, gdcm\nimport torch\nfrom torch import nn\nfrom torch.nn import functional as F\nfrom torchvision.models import resnet18, ResNet18_Weights\nfrom sklearn.model_selection import KFold\nfrom sklearn.metrics import roc_auc_score, roc_curve, average_precision_score\nfrom tqdm.auto import tqdm\nfrom IPython.display import display, Markdown\n\nDEVICE = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\ndef seed_all(seed):\n    random.seed(seed); np.random.seed(seed); torch.manual_seed(seed)\n    if torch.cuda.is_available(): torch.cuda.manual_seed_all(seed)\n    torch.backends.cudnn.benchmark = False\n    torch.backends.cudnn.deterministic = True\nseed_all(CFG[\"seed\"])\nversions = {p: metadata.version(p) for p in [\"torch\", \"torchvision\", \"numpy\", \"pandas\", \"scikit-learn\", \"pydicom\", \"python-gdcm\"]}\ndisplay(pd.Series(versions, name=\"version\").to_frame())\nprint(\"裝置：\", DEVICE, \"GPU 數量：\", torch.cuda.device_count())\nif DEVICE.type == \"cpu\": print(\"目前沒有 GPU；影像特徵抽取會比較慢。\")\nif MODE == \"train\":\n    (OUT / \"wheels\").mkdir(exist_ok=True)\n    for dist in packages.values():\n        subprocess.check_call([sys.executable, \"-m\", \"pip\", \"download\", \"--quiet\", \"--no-deps\",\n                               \"--only-binary=:all:\", \"-d\", str(OUT / \"wheels\"), f\"{dist}=={metadata.version(dist)}\"])\n    (OUT / \"environment.json\").write_text(json.dumps(versions, indent=2))\n","id":"knee-005"},{"cell_type":"markdown","metadata":{},"source":"## 04｜確認資料，而不是直接開始訓練\n\n我先檢查欄位、重複 study 與標籤值。若資料更新，程式會明確停止，不會偷偷選其他欄位。缺失與 -1 先保留為未知；其他非 0/1 數值必須回頭核對。\n","id":"knee-006"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"key = \"train.csv\" if MODE == \"train\" else \"test.csv\"\nroots = [p.parent for p in INPUT.rglob(key)\n         if (p.parent / (\"train_series.csv\" if MODE == \"train\" else \"test_series.csv\")).exists()\n         and (p.parent / (\"train_series\" if MODE == \"train\" else \"test_series\")).is_dir()]\nif len(roots) != 1:\n    raise RuntimeError(f\"找到 {len(roots)} 個比賽目錄；請確認只掛載一份官方比賽資料。\")\nROOT = roots[0]\nsplit_name = \"train\" if MODE == \"train\" else \"test\"\nstudies = pd.read_csv(ROOT / f\"{split_name}.csv\", dtype={\"StudyInstanceUID\": str})\nseries = pd.read_csv(ROOT / f\"{split_name}_series.csv\", dtype={\"StudyInstanceUID\": str, \"SeriesInstanceUID\": str})\nrequired = {\"StudyInstanceUID\", \"SeriesInstanceUID\", \"Fluid_Sensitive\", \"Fat_Suppression\", \"Anatomical_Plane\"}\nassert required <= set(series.columns), f\"序列欄位不同：{series.columns.tolist()}\"\nassert studies[\"StudyInstanceUID\"].notna().all() and not studies[\"StudyInstanceUID\"].duplicated().any()\nassert not series.duplicated([\"StudyInstanceUID\", \"SeriesInstanceUID\"]).any()\nprint(\"資料目錄：\", ROOT, \"｜檢查數：\", len(studies), \"｜序列數：\", len(series))\ndisplay(studies.drop(columns=[\"Report\"], errors=\"ignore\").head())\ndisplay(series.head())\nif MODE == \"train\":\n    assert set(TARGETS) <= set(studies.columns)\n    for c in TARGETS:\n        studies[c] = pd.to_numeric(studies[c], errors=\"raise\").replace(-1, np.nan)\n        assert studies[c].dropna().isin([0, 1]).all(), f\"{c} 出現非二元標籤\"\n    gold = studies[studies[TARGETS].notna().any(axis=1)].reset_index(drop=True)\n    assert len(gold) >= CFG[\"folds\"] * 2, \"人工標籤樣本不足以進行目前的交叉驗證。\"\n    print(\"至少一個人工標籤：\", len(gold), \"｜完全未標註：\", len(studies)-len(gold))\n    print(\"人工標籤完整的檢查：\", studies[TARGETS].notna().all(axis=1).sum())\n","id":"knee-007"},{"cell_type":"markdown","metadata":{},"source":"## 05｜標籤數量與資料分布\n\n我同時看陽性、陰性與缺失，避免把樣本數差異藏起來。下圖只統計人工標籤，不從報告猜答案。圖內使用官方英文欄位，避免 Kaggle 缺少中文字型導致方框。\n","id":"knee-008"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"FIG = OUT / \"figures\"\nFIG.mkdir(exist_ok=True)\ndef savefig(name):\n    plt.tight_layout(); plt.savefig(FIG / name, dpi=150, bbox_inches=\"tight\"); plt.show()\nif MODE == \"train\":\n    counts = pd.DataFrame({\"positive\": (gold[TARGETS] == 1).sum(),\n                           \"negative\": (gold[TARGETS] == 0).sum(),\n                           \"unknown\": gold[TARGETS].isna().sum()})\n    display(counts)\n    counts.to_csv(OUT / \"label_counts.csv\")\n    counts.plot.barh(stacked=True, figsize=(10, 6), color=[\"#e49b53\", \"#4b86a8\", \"#b9bfc5\"])\n    plt.title(\"Reviewed labels only\"); plt.xlabel(\"Study-label count\")\n    savefig(\"01_label_distribution.png\")\n","id":"knee-009"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    fig, ax = plt.subplots(1, 2, figsize=(11, 4))\n    series.groupby(\"StudyInstanceUID\").size().plot.hist(bins=20, ax=ax[0])\n    ax[0].set(title=\"Series per study\", xlabel=\"Number of series\")\n    series[\"Anatomical_Plane\"].fillna(\"Unknown\").value_counts().plot.bar(ax=ax[1])\n    ax[1].set(title=\"Anatomical planes\", ylabel=\"Series count\")\n    savefig(\"02_series_distribution.png\")\n    if \"Report\" in studies:\n        lengths = studies[\"Report\"].fillna(\"\").str.len()\n        display(lengths.describe().to_frame(\"Report characters\"))\n        lengths.plot.hist(bins=40, figsize=(8, 3)); plt.title(\"Report length (characters)\")\n        savefig(\"03_report_lengths.png\")\n        print(\"這裡只量測報告長度，不把文字長短或語言當成病況標籤。\")\n","id":"knee-010"},{"cell_type":"markdown","metadata":{},"source":"## 06｜我如何挑選 MRI 序列\n\n每次檢查可能包含很多序列。我先在矢狀面、冠狀面、軸位各挑一個，優先選 fluid-sensitive、其次 fat-suppressed 的序列，並用 UID 排序確保結果可重現。這是基準版的選片規則，不代表它對每種病灶都最佳。\n\n切片依 ImageOrientationPatient 與 ImagePositionPatient 的幾何投影排序；缺少幾何資訊才改用 InstanceNumber。每個切面均勻取 8 張，最多 24 張。若不足 8 張，只保留不同切片並遮罩空位，不重複假裝有更多影像。這版未做完整的物理空間重採樣或左右翻轉標準化，跨院方向差異仍是限制。\n","id":"knee-011"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"def choose_series(frame):\n    frame = frame.copy()\n    for c in [\"Fluid_Sensitive\", \"Fat_Suppression\"]:\n        frame[c] = pd.to_numeric(frame[c], errors=\"coerce\").fillna(0)\n    frame[\"Anatomical_Plane\"] = frame[\"Anatomical_Plane\"].astype(str).str.strip().str.title()\n    frame = frame[frame[\"Anatomical_Plane\"].isin(PLANES)]\n    frame = frame.sort_values([\"StudyInstanceUID\", \"Anatomical_Plane\", \"Fluid_Sensitive\",\n                               \"Fat_Suppression\", \"SeriesInstanceUID\"], ascending=[True, True, False, False, True])\n    return frame.drop_duplicates([\"StudyInstanceUID\", \"Anatomical_Plane\"])\nchosen = choose_series(series)\nseries_map = {uid: g for uid, g in chosen.groupby(\"StudyInstanceUID\")}\ndisplay(chosen.head(9))\nprint(\"各檢查可用切面數：\")\ndisplay(chosen.groupby(\"StudyInstanceUID\").size().value_counts().sort_index())\n\ndef sorted_dicoms(folder):\n    paths = sorted(folder.glob(\"*.dcm\"))\n    if not paths: raise RuntimeError(f\"沒有 DICOM：{folder}\")\n    headers = [pydicom.dcmread(p, stop_before_pixels=True) for p in paths]\n    geometry_ok = all(hasattr(h, \"ImageOrientationPatient\") and hasattr(h, \"ImagePositionPatient\") for h in headers)\n    if geometry_ok:\n        ori = np.array(headers[0].ImageOrientationPatient, dtype=float)\n        normal = np.cross(ori[:3], ori[3:])\n        geometry_ok = np.linalg.norm(normal) > 0.9 and all(\n            np.allclose(np.array(h.ImageOrientationPatient, float), ori, atol=0.02) for h in headers)\n    if geometry_ok:\n        keys = [float(np.dot(np.array(h.ImagePositionPatient, float), normal)) for h in headers]\n    elif all(hasattr(h, \"InstanceNumber\") for h in headers):\n        keys = [float(h.InstanceNumber) for h in headers]\n    else:\n        raise RuntimeError(f\"無法可靠排序：{folder}；需要檢查 DICOM metadata。\")\n    order = sorted(range(len(paths)), key=lambda i: (keys[i], str(paths[i])))\n    return [paths[i] for i in order], headers[order[0]], \"geometry\" if geometry_ok else \"instance\"\n\ndef read_slice(path):\n    ds = pydicom.dcmread(path)\n    try:\n        x = ds.pixel_array.astype(np.float32)\n    except Exception as e:\n        raise RuntimeError(f\"DICOM 解碼失敗：{path.name}，syntax={ds.file_meta.TransferSyntaxUID}；檢查 GDCM。\") from e\n    if x.ndim != 2: raise RuntimeError(f\"預期單張灰階影像，實際 {x.shape}\")\n    x = x * float(getattr(ds, \"RescaleSlope\", 1)) + float(getattr(ds, \"RescaleIntercept\", 0))\n    if not np.isfinite(x).all(): raise RuntimeError(f\"像素包含非有限值：{path.name}\")\n    # 每張使用 1–99 百分位裁切；不使用整份資料的統計，避免驗證資料洩漏。\n    lo, hi = np.percentile(x, [1, 99])\n    x = np.clip((x-lo) / max(float(hi-lo), 1e-6), 0, 1)\n    if getattr(ds, \"PhotometricInterpretation\", \"\") == \"MONOCHROME1\": x = 1-x\n    h, w = x.shape\n    scale = CFG[\"image_size\"] / max(h, w)\n    nh, nw = max(1, round(h*scale)), max(1, round(w*scale))\n    t = F.interpolate(torch.from_numpy(x)[None, None], size=(nh, nw), mode=\"bilinear\", align_corners=False)[0]\n    dh, dw = CFG[\"image_size\"]-nh, CFG[\"image_size\"]-nw\n    return F.pad(t, (dw//2, dw-dw//2, dh//2, dh-dh//2))\n\ndef load_study(uid, which):\n    if uid not in series_map: raise RuntimeError(f\"檢查 {uid} 沒有三個指定切面的任何序列。\")\n    n = CFG[\"slices_per_plane\"]\n    images = torch.zeros(3*n, 1, CFG[\"image_size\"], CFG[\"image_size\"])\n    mask = np.zeros(3*n, dtype=bool)\n    info, patient_ids, scanners = [], set(), set()\n    for row in series_map[uid].itertuples(index=False):\n        plane = row.Anatomical_Plane; pi = PLANES.index(plane)\n        folder = ROOT / f\"{which}_series\" / uid / row.SeriesInstanceUID\n        paths, header, method = sorted_dicoms(folder)\n        pid = str(getattr(header, \"PatientID\", \"\")).strip()\n        issuer = str(getattr(header, \"IssuerOfPatientID\", \"\")).strip()\n        if pid: patient_ids.add(issuer + \"::\" + pid)\n        scanners.add(str(getattr(header, \"Manufacturer\", \"unknown\")))\n        selected = np.unique(np.linspace(0, len(paths)-1, min(n, len(paths))).round().astype(int))\n        for j, index in enumerate(selected):\n            slot = pi*n+j\n            images[slot] = read_slice(paths[index]); mask[slot] = True\n            info.append(dict(slot=slot, plane=plane, slice_index=int(index), total_slices=len(paths),\n                             sort_method=method, file=paths[index].name))\n    if len(patient_ids) > 1: raise RuntimeError(f\"同一 study 有不同 PatientID：{uid}\")\n    # 雜湊只做分組，不輸出原始 PatientID。\n    group = \"patient:\" + hashlib.sha256(next(iter(patient_ids)).encode()).hexdigest() if patient_ids else \"study:\"+uid\n    return images, mask, info, group, \"|\".join(sorted(scanners))\n","id":"knee-012"},{"cell_type":"markdown","metadata":{},"source":"## 07｜先看一個真實樣本\n\n這裡的影像由資料實際讀出，不是示意圖。我先確認三個切面、亮暗與切片順序是否合理。空白格代表沒有該切片，不會參與模型聚合。\n","id":"knee-013"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    example_uid = gold.iloc[0][\"StudyInstanceUID\"]\n    example_x, example_mask, example_info, _, _ = load_study(example_uid, \"train\")\n    fig, axes = plt.subplots(3, CFG[\"slices_per_plane\"], figsize=(16, 6), squeeze=False)\n    for i, ax in enumerate(axes.flat):\n        ax.imshow(example_x[i, 0], cmap=\"gray\", vmin=0, vmax=1); ax.axis(\"off\")\n        ax.set_title(f\"{PLANES[i//CFG['slices_per_plane']]} {i%CFG['slices_per_plane']+1}\" if example_mask[i] else \"Missing\", fontsize=8)\n    savefig(\"04_mri_contact_sheet.png\")\n    display(pd.DataFrame(example_info).drop(columns=\"file\"))\n    display(gold.iloc[0][TARGETS].to_frame(\"Reviewed label\"))\n","id":"knee-014"},{"cell_type":"markdown","metadata":{},"source":"## 08｜我先用預訓練 CNN 擷取特徵\n\n我使用 ImageNet 預訓練 ResNet-18，把每張切片轉成 512 維特徵；灰階複製為三通道，再套用 ImageNet normalization。為保留完整膝蓋，我採等比例縮放加補邊，並未使用原始 ImageNet 的中心裁切。\n\n第一版固定 CNN，不微調它，只訓練後面的聚合與分類模組。這樣能減少小樣本過擬合與 GPU 負擔，但 ImageNet 與 MRI 的領域差距仍存在。特徵可以一次快取，兩個實驗共用；CNN 沒有用到本比賽標籤，因此在切分前抽取固定特徵不會讓它學到驗證答案。\n","id":"knee-015"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"def atomic_torch_save(obj, path):\n    path = Path(path); tmp = path.with_suffix(path.suffix + \".tmp\")\n    torch.save(obj, tmp); os.replace(tmp, path)\n\nbackbone = resnet18(weights=None)\nif MODE == \"infer\":\n    backbone.fc = nn.Identity()\n    backbone.load_state_dict(torch.load(SOURCE / \"backbone.pt\", map_location=\"cpu\", weights_only=True))\nelif SOURCE is not None and (SOURCE / \"backbone.pt\").exists():\n    backbone.fc = nn.Identity()\n    backbone.load_state_dict(torch.load(SOURCE / \"backbone.pt\", map_location=\"cpu\", weights_only=True))\nelse:\n    # 第一次訓練需 Internet；失敗時停止，不偷偷改成隨機權重。\n    backbone = resnet18(weights=ResNet18_Weights.IMAGENET1K_V1)\n    backbone.fc = nn.Identity()\nbackbone = backbone.to(DEVICE).eval()\nfor p in backbone.parameters(): p.requires_grad_(False)\nif MODE == \"train\": atomic_torch_save(backbone.cpu().state_dict(), OUT / \"backbone.pt\"); backbone.to(DEVICE)\n\ndef state_hash(model):\n    h = hashlib.sha256()\n    for k, v in model.state_dict().items(): h.update(k.encode()); h.update(v.detach().cpu().numpy().tobytes())\n    return h.hexdigest()\nBACKBONE_HASH = state_hash(backbone)\nFEATURE_SIGNATURE = hashlib.sha256(json.dumps([CFG[\"model_version\"], CFG[\"image_size\"],\n                                              CFG[\"slices_per_plane\"], BACKBONE_HASH]).encode()).hexdigest()\n@torch.inference_mode()\ndef extract_features(images, mask):\n    features = torch.zeros(len(mask), 512)\n    valid = np.flatnonzero(mask)\n    for start in range(0, len(valid), CFG[\"feature_batch_size\"]):\n        ids = valid[start:start+CFG[\"feature_batch_size\"]]\n        x = images[ids].repeat(1, 3, 1, 1).to(DEVICE)\n        mean = x.new_tensor([0.485, 0.456, 0.406])[None, :, None, None]\n        std = x.new_tensor([0.229, 0.224, 0.225])[None, :, None, None]\n        features[ids] = backbone((x-mean)/std).float().cpu()\n    return features.numpy()\n\nCACHE = OUT / \"features\"\nCACHE.mkdir(exist_ok=True)\ndef get_features(uid, which):\n    filename = f\"{which}_{uid}.npz\"\n    candidates = [CACHE / filename]\n    if SOURCE is not None: candidates.append(SOURCE / \"features\" / filename)\n    for path in candidates:\n        if path.exists():\n            with np.load(path, allow_pickle=False) as z:\n                if str(z[\"signature\"]) == FEATURE_SIGNATURE:\n                    result = tuple(z[k].copy() for k in [\"features\", \"mask\", \"group\", \"scanner\"])\n                    if path != CACHE / filename: shutil.copy2(path, CACHE / filename)\n                    return result\n    x, mask, info, group, scanner = load_study(uid, which)\n    features = extract_features(x, mask)\n    path = CACHE / filename\n    with open(str(path)+\".tmp\", \"wb\") as f:\n        np.savez_compressed(f, features=features, mask=mask, group=group, scanner=scanner, signature=FEATURE_SIGNATURE)\n    os.replace(str(path)+\".tmp\", path)\n    return features, mask, np.array(group), np.array(scanner)\n","id":"knee-016"},{"cell_type":"markdown","metadata":{},"source":"## 09｜建立可續用的特徵快取\n\n每完成一個 study 就保存一個小型特徵檔。遇到壞檔會停止並顯示位置，不用全黑影像或 0.5 偷偷補掉。第一次跑完後，我會用實測耗時估計正式測試成本；有些壓縮格式或長序列可能比較慢，估計值不是保證。\n","id":"knee-017"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    start = time.monotonic()\n    extracted = [get_features(uid, \"train\") for uid in tqdm(gold.StudyInstanceUID, desc=\"人工標籤影像特徵\")]\n    X = np.stack([x[0] for x in extracted]).astype(np.float32)\n    M = np.stack([x[1] for x in extracted]).astype(bool)\n    GROUPS = np.array([str(x[2]) for x in extracted])\n    SCANNERS = np.array([str(x[3]) for x in extracted])\n    Y = gold[TARGETS].to_numpy(dtype=np.float32)\n    print(\"特徵張量：\", X.shape, \"｜有效切片比例：\", M.mean().round(3))\n    print(\"本次抽取／讀取快取耗時（秒）：\", round(time.monotonic()-start, 1))\n    print(\"病人識別可用的檢查數：\", sum(g.startswith(\"patient:\") for g in GROUPS))\n    assert np.isfinite(X).all() and M.any(axis=1).all()\n","id":"knee-018"},{"cell_type":"markdown","metadata":{},"source":"## 10｜切分的單位是檢查／病人，不是切片\n\n同一次檢查的切片不能同時出現在訓練與驗證。DICOM 若提供 PatientID，我把相同 ID 放在同一 fold；沒有時只能以 study 分組。匿名 ID 未必能完整連結跨檢查的同一病人，因此我不能宣稱已完全排除病人重複。這版也沒有把醫院完整隔離，不能把分數當成跨院驗證。\n\n我先固定 5 folds；不為了分數反覆挑 seed。小樣本下某 fold 的某類可能只有陽性或陰性，此時 AUC 不存在，我記為 NaN。最後集中所有 out-of-fold（OOF）預測計算每類 AUC，並顯示實際可計算的類別數。\n","id":"knee-019"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    unique_groups = np.unique(GROUPS)\n    assert len(unique_groups) >= CFG[\"folds\"], \"可用病人／study 分組不足；需檢查匿名 ID，不能直接拆散。\"\n    fold_id = np.full(len(gold), -1, dtype=int)\n    splitter = KFold(n_splits=CFG[\"folds\"], shuffle=True, random_state=CFG[\"seed\"])\n    for fold, (_, vi) in enumerate(splitter.split(unique_groups)):\n        fold_id[np.isin(GROUPS, unique_groups[vi])] = fold\n    assert (fold_id >= 0).all()\n    for fold in range(CFG[\"folds\"]):\n        assert not (set(GROUPS[fold_id==fold]) & set(GROUPS[fold_id!=fold]))\n    split_table = gold[[\"StudyInstanceUID\"]].assign(group=GROUPS, manufacturer=SCANNERS, fold=fold_id)\n    split_table.to_csv(OUT / \"folds.csv\", index=False)\n    display(pd.crosstab(split_table.fold, split_table.manufacturer))\n    display(pd.DataFrame([{ \"fold\": f, \"studies\": int((fold_id==f).sum()),\n                            **{c: int(np.nansum(Y[fold_id==f, j])) for j, c in enumerate(TARGETS)}}\n                         for f in range(CFG[\"folds\"])]))\n","id":"knee-020"},{"cell_type":"markdown","metadata":{},"source":"## 11｜我的對照實驗\n\n| 實驗 | 相同部分 | 唯一主要差異 |\n|---|---|---|\n| mean | 固定 CNN、相同切片、相同 folds、相同訓練輪數 | 有效切片等權平均 |\n| attention | 同上 | 以可學習權重聚合切片 |\n\n兩者都加入切面 embedding，避免模型完全不知道切片來自哪個方向。注意力權重是對本模型決策的相對加權，不是病灶位置標註，也不是醫師認可的解釋。這裡使用共用注意力，不是每種疾病各自一套注意力。\n\n我事先指定 `attention` 用於提交，不根據這份小樣本 OOF 的高低自動挑選模型。比較結果仍屬探索性分析，未另設獨立測試集。\n","id":"knee-021"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"class KneeMIL(nn.Module):\n    def __init__(self, pooling):\n        super().__init__()\n        assert pooling in {\"mean\", \"attention\"}\n        self.pooling = pooling\n        self.proj = nn.Sequential(nn.LayerNorm(512), nn.Linear(512, 128), nn.GELU())\n        self.plane_embedding = nn.Embedding(3, 128)\n        self.attention = nn.Sequential(nn.Linear(128, 64), nn.Tanh(), nn.Linear(64, 1))\n        self.head = nn.Sequential(nn.Dropout(CFG[\"dropout\"]), nn.Linear(128, len(TARGETS)))\n    def forward(self, x, mask, return_attention=False):\n        plane = torch.arange(3, device=x.device).repeat_interleave(CFG[\"slices_per_plane\"])\n        h = self.proj(x) + self.plane_embedding(plane)[None]\n        if self.pooling == \"attention\":\n            scores = self.attention(h).squeeze(-1).masked_fill(~mask, -1e4)\n            weights = torch.softmax(scores, dim=1)\n        else:\n            weights = mask.float()/mask.sum(1, keepdim=True).clamp_min(1)\n        logits = self.head((h*weights[..., None]).sum(1))\n        return (logits, weights) if return_attention else logits\n\ndef masked_bce(logits, target):\n    known = torch.isfinite(target)\n    loss = F.binary_cross_entropy_with_logits(logits, torch.nan_to_num(target, nan=0.0), reduction=\"none\")\n    return (loss*known).sum()/known.sum().clamp_min(1)\n\ndef metric_table(y, p):\n    result = []\n    for j, c in enumerate(TARGETS):\n        ok = np.isfinite(y[:, j]); yt, pt = y[ok, j], p[ok, j]\n        valid = len(np.unique(yt)) == 2\n        result.append(dict(target=c, n=len(yt), positives=int(yt.sum()),\n                           auc=roc_auc_score(yt, pt) if valid else np.nan,\n                           ap=average_precision_score(yt, pt) if valid else np.nan,\n                           brier=float(np.mean((yt-pt)**2)) if len(yt) else np.nan))\n    return pd.DataFrame(result)\n\n@torch.inference_mode()\ndef predict_head(model, x, mask, batch=32):\n    model.eval(); output=[]\n    for start in range(0, len(x), batch):\n        xb = torch.as_tensor(x[start:start+batch], device=DEVICE)\n        mb = torch.as_tensor(mask[start:start+batch], device=DEVICE)\n        output.append(model(xb, mb).sigmoid().cpu().numpy())\n    return np.concatenate(output)\n","id":"knee-022"},{"cell_type":"markdown","metadata":{},"source":"## 12｜訓練、紀錄與 checkpoint\n\n我固定訓練 40 epochs，保存最後一輪模型；不根據驗證表現挑最佳 epoch，因此 OOF 不包含挑 epoch 的額外偏差。每一輪記錄 loss 與可計算類別的 AUC。訓練只對已知標籤計算 BCE，未知標籤不計入 loss。\n\ncheckpoint 包含 optimizer 與亂數狀態，且檢查資料／設定簽章。載入使用 `weights_only=True`，只讀取本實驗格式。執行速度主要取決於前面的 DICOM 讀取與 CNN，這裡的分類頭通常較快。\n","id":"knee-023"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"def experiment_signature():\n    h=hashlib.sha256()\n    h.update(json.dumps(CFG, sort_keys=True).encode()); h.update(FEATURE_SIGNATURE.encode())\n    h.update(gold.StudyInstanceUID.str.cat(sep=\"|\").encode())\n    h.update(Y.tobytes()); h.update(fold_id.tobytes()); h.update(X.tobytes()); h.update(M.tobytes())\n    return h.hexdigest()\n\ndef train_fold(pooling, fold, signature):\n    seed_all(CFG[\"seed\"]+fold)\n    model = KneeMIL(pooling).to(DEVICE)\n    optimizer = torch.optim.AdamW(model.parameters(), lr=CFG[\"lr\"], weight_decay=CFG[\"weight_decay\"])\n    tr, va = np.flatnonzero(fold_id != fold), np.flatnonzero(fold_id == fold)\n    path = OUT / f\"{pooling}_fold{fold}.pt\"\n    old = path if path.exists() else (SOURCE / path.name if SOURCE is not None else path)\n    begin, logs = 0, []\n    if old.exists():\n        ck = torch.load(old, map_location=\"cpu\", weights_only=True)\n        if ck[\"signature\"] != signature: raise RuntimeError(\"舊 checkpoint 設定／資料不同；請改用新的 OUT。\")\n        model.load_state_dict(ck[\"state\"]); optimizer.load_state_dict(ck[\"optimizer\"])\n        for state in optimizer.state.values():\n            for k, v in state.items():\n                if torch.is_tensor(v): state[k] = v.to(DEVICE)\n        begin, logs = ck[\"epoch\"], ck[\"logs\"]\n        torch.set_rng_state(ck[\"rng\"])\n        if DEVICE.type == \"cuda\" and ck[\"cuda_rng\"]: torch.cuda.set_rng_state_all(ck[\"cuda_rng\"])\n        if old != path: shutil.copy2(old, path)\n    xt, mt, yt = (torch.as_tensor(a[tr], device=DEVICE) for a in [X, M, Y])\n    for epoch in range(begin, CFG[\"epochs\"]):\n        model.train(); total, batches = 0., 0\n        order = torch.randperm(len(tr), device=DEVICE)\n        for start in range(0, len(tr), CFG[\"batch_size\"]):\n            ids = order[start:start+CFG[\"batch_size\"]]\n            optimizer.zero_grad(set_to_none=True)\n            loss = masked_bce(model(xt[ids], mt[ids]), yt[ids])\n            if not torch.isfinite(loss): raise RuntimeError(\"Loss 非有限值；停止保存。\")\n            loss.backward(); nn.utils.clip_grad_norm_(model.parameters(), 1.0); optimizer.step()\n            total += loss.item(); batches += 1\n        pred = predict_head(model, X[va], M[va])\n        metrics = metric_table(Y[va], pred)\n        known = np.isfinite(Y[va]); py = np.clip(pred, 1e-7, 1-1e-7)\n        vy = np.nan_to_num(Y[va], nan=0.)\n        val_loss = float((-(vy*np.log(py)+(1-vy)*np.log(1-py)))[known].mean())\n        row = dict(pooling=pooling, fold=fold, epoch=epoch+1, train_loss=total/max(batches,1),\n                   val_loss=val_loss, val_auc=float(metrics.auc.mean()), valid_targets=int(metrics.auc.notna().sum()))\n        logs.append(row)\n        atomic_torch_save(dict(state={k:v.detach().cpu() for k,v in model.state_dict().items()},\n                               optimizer=optimizer.state_dict(), epoch=epoch+1, logs=logs,\n                               signature=signature, pooling=pooling, rng=torch.get_rng_state(),\n                               cuda_rng=torch.cuda.get_rng_state_all() if DEVICE.type==\"cuda\" else []), path)\n        if epoch == begin or (epoch+1)%10 == 0: print(row, flush=True)\n    pred = predict_head(model, X[va], M[va])\n    del model, optimizer, xt, mt, yt\n    if DEVICE.type==\"cuda\": torch.cuda.empty_cache()\n    return va, pred, logs\n","id":"knee-024"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    SIGNATURE = experiment_signature()\n    OOF, history_rows = {}, []\n    for pooling in [\"mean\", \"attention\"]:\n        oof = np.full_like(Y, np.nan)\n        for fold in range(CFG[\"folds\"]):\n            va, pred, logs = train_fold(pooling, fold, SIGNATURE)\n            oof[va] = pred; history_rows.extend(logs)\n        assert np.isfinite(oof).all()\n        OOF[pooling] = oof\n        frame = pd.DataFrame(oof, columns=TARGETS)\n        frame.insert(0, \"StudyInstanceUID\", gold.StudyInstanceUID)\n        frame.to_csv(OUT / f\"oof_{pooling}.csv\", index=False)\n    history = pd.DataFrame(history_rows)\n    history.to_csv(OUT / \"history.csv\", index=False)\n    print(\"兩組交叉驗證完成；以下所有驗證預測都來自未看過該分組的模型。\")\n","id":"knee-025"},{"cell_type":"markdown","metadata":{},"source":"## 13｜我先看訓練曲線\n\n如果 train loss 持續下降，validation loss 卻上升，我會記錄過擬合，而不是只留下看起來最好的一張圖。各 fold 的 validation AUC 可能劇烈波動；在樣本少的情況下，這不一定代表模型穩定進步。\n","id":"knee-026"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    fig, axes = plt.subplots(1, 2, figsize=(12, 4))\n    for pooling, g in history.groupby(\"pooling\"):\n        avg = g.groupby(\"epoch\")[[\"train_loss\", \"val_loss\", \"val_auc\"]].mean()\n        axes[0].plot(avg.index, avg.train_loss, label=pooling+\" train\")\n        axes[0].plot(avg.index, avg.val_loss, \"--\", label=pooling+\" validation\")\n        axes[1].plot(avg.index, avg.val_auc, label=pooling)\n    axes[0].set(title=\"Mean loss across folds\", xlabel=\"Epoch\", ylabel=\"BCE\")\n    axes[1].set(title=\"Mean fold AUC (valid targets)\", xlabel=\"Epoch\", ylabel=\"ROC-AUC\")\n    for ax in axes: ax.legend()\n    savefig(\"05_learning_curves.png\")\n","id":"knee-027"},{"cell_type":"markdown","metadata":{},"source":"## 14｜OOF 評估與簡單基線\n\n我另外加入「每類都預測訓練 fold 盛行率」的無影像基線。跨 fold 的盛行率不同，合併 OOF 的常數基線 AUC 不一定剛好是 0.5，因此我也保留每類結果。當有類別完全無法計算 AUC 時，顯示的平均只涵蓋有效類別，不能直接當成官方完整 12 類分數。\n\n除了 AUC，我也看 AP（precision-recall 的摘要）與 Brier score（機率誤差）。sigmoid 輸出只是未校準的模型分數，不能直接解讀成實際患病機率。\n","id":"knee-028"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    prior_oof = np.zeros_like(Y)\n    for f in range(CFG[\"folds\"]):\n        yt = Y[fold_id!=f]\n        # Laplace smoothing，分母只計已知標籤。\n        prior_oof[fold_id==f] = (np.nansum(yt, axis=0)+1)/(np.isfinite(yt).sum(axis=0)+2)\n    OOF[\"prior\"] = prior_oof\n    all_metrics=[]\n    for name, p in OOF.items():\n        table=metric_table(Y, p); table[\"experiment\"]=name; all_metrics.append(table)\n    metrics=pd.concat(all_metrics, ignore_index=True)\n    metrics.to_csv(OUT / \"metrics_by_target.csv\", index=False)\n    summary=metrics.groupby(\"experiment\").agg(macro_auc=(\"auc\",\"mean\"), valid_targets=(\"auc\",\"count\"),\n                                              macro_ap=(\"ap\",\"mean\"), mean_brier=(\"brier\",\"mean\"))\n    summary.to_csv(OUT / \"experiment_summary.csv\")\n    display(summary)\n    display(metrics.pivot(index=\"target\", columns=\"experiment\", values=\"auc\"))\n    metrics.pivot(index=\"target\", columns=\"experiment\", values=\"auc\").plot.barh(figsize=(10,7))\n    plt.axvline(0.5, color=\"gray\", linestyle=\"--\"); plt.xlim(0,1); plt.xlabel(\"OOF ROC-AUC\")\n    savefig(\"06_auc_comparison.png\")\n","id":"knee-029"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    fig, axes=plt.subplots(3,4, figsize=(14,10))\n    for j, (name, ax) in enumerate(zip(TARGETS, axes.flat)):\n        ok=np.isfinite(Y[:,j])\n        for pooling in [\"mean\",\"attention\"]:\n            if np.unique(Y[ok,j]).size==2:\n                fpr,tpr,_=roc_curve(Y[ok,j], OOF[pooling][ok,j]); ax.plot(fpr,tpr,label=pooling)\n        ax.plot([0,1],[0,1],\"--\",color=\"gray\"); ax.set(title=name, xlabel=\"False positive rate\", ylabel=\"True positive rate\")\n    axes[0,0].legend(); savefig(\"07_roc_curves.png\")\n","id":"knee-030"},{"cell_type":"markdown","metadata":{},"source":"## 15｜小樣本分數有多不穩定？\n\n我以病人／study 分組重抽樣，對同一組 OOF 預測估計 bootstrap 區間，同時比較 attention − mean 的配對差異。每次抽樣要求原本可評估的類別都仍有正反例；不滿足就跳過。這個區間只反映有限樣本在既有 OOF 預測上的波動，**沒有涵蓋重新訓練、換資料來源或模型選擇的變異**。\n","id":"knee-031"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"def bootstrap_oof(y, a, b, groups, repeats, seed):\n    rng=np.random.default_rng(seed); unique=np.unique(groups)\n    valid=[j for j in range(y.shape[1]) if np.unique(y[np.isfinite(y[:,j]),j]).size==2]\n    values=[]\n    if not valid: return pd.DataFrame(columns=[\"attention_auc\",\"attention_minus_mean\"])\n    blocks={g:np.flatnonzero(groups==g) for g in unique}\n    for _ in range(repeats):\n        ids=np.concatenate([blocks[g] for g in rng.choice(unique,len(unique),replace=True)])\n        aa,bb=[],[]\n        for j in valid:\n            ok=np.isfinite(y[ids,j]); ii=ids[ok]\n            if np.unique(y[ii,j]).size!=2: break\n            aa.append(roc_auc_score(y[ii,j],a[ii,j])); bb.append(roc_auc_score(y[ii,j],b[ii,j]))\n        if len(aa)==len(valid): values.append((np.mean(aa),np.mean(aa)-np.mean(bb)))\n    return pd.DataFrame(values,columns=[\"attention_auc\",\"attention_minus_mean\"])\n\nif MODE == \"train\":\n    boot=bootstrap_oof(Y, OOF[\"attention\"], OOF[\"mean\"], GROUPS, CFG[\"bootstrap\"], CFG[\"seed\"])\n    boot.to_csv(OUT / \"bootstrap.csv\", index=False)\n    print(\"有效 bootstrap 次數：\", len(boot), \"/\", CFG[\"bootstrap\"])\n    if len(boot)>=30:\n        display(boot.quantile([0.025,0.5,0.975]))\n        boot.attention_minus_mean.plot.hist(bins=25,figsize=(7,3))\n        plt.axvline(0,color=\"red\",linestyle=\"--\"); plt.title(\"Paired bootstrap: attention minus mean AUC\")\n        savefig(\"08_bootstrap_difference.png\")\n    else: print(\"有效重抽樣太少，我暫時不把區間當成穩定估計。\")\n","id":"knee-032"},{"cell_type":"markdown","metadata":{},"source":"## 16｜錯誤案例與切片權重\n\n我找出 OOF 平均機率誤差最大的檢查，再載入當時沒有看過它的 fold 模型，避免拿訓練樣本做漂亮的展示。紅／藍色熱圖是「標籤與模型分數」，不是 MRI 上的病灶分割。\n","id":"knee-033"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    errors=np.nanmean(np.abs(OOF[\"attention\"]-Y),axis=1)\n    worst=np.argsort(-errors)[:min(5,len(Y))]\n    error_table=gold.iloc[worst][[\"StudyInstanceUID\"]].copy()\n    error_table[\"mean_absolute_probability_error\"]=errors[worst]\n    error_table.to_csv(OUT / \"error_cases.csv\",index=False); display(error_table)\n    case=int(worst[0]); f=int(fold_id[case])\n    review_model=KneeMIL(\"attention\").to(DEVICE)\n    review_ck=torch.load(OUT / f\"attention_fold{f}.pt\",map_location=\"cpu\",weights_only=True)\n    review_model.load_state_dict(review_ck[\"state\"]); review_model.eval()\n    with torch.inference_mode():\n        logits, weights=review_model(torch.as_tensor(X[case:case+1],device=DEVICE),\n                                    torch.as_tensor(M[case:case+1],device=DEVICE),True)\n    np.testing.assert_allclose(logits.sigmoid().cpu().numpy()[0], OOF[\"attention\"][case], atol=1e-5)\n    plt.figure(figsize=(12,2))\n    plt.imshow(np.vstack([Y[case], OOF[\"attention\"][case]]), vmin=0,vmax=1,cmap=\"coolwarm\",aspect=\"auto\")\n    plt.xticks(range(12),TARGETS,rotation=45,ha=\"right\"); plt.yticks([0,1],[\"Reviewed label\",\"OOF score\"])\n    plt.colorbar(); savefig(\"09_error_case_scores.png\")\n","id":"knee-034"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    cx, cm, ci, _, _=load_study(gold.iloc[case].StudyInstanceUID,\"train\")\n    w=weights[0].cpu().numpy()\n    top=np.argsort(-np.where(cm,w,-1))[:min(6,cm.sum())]\n    fig, axes=plt.subplots(1,len(top),figsize=(15,3),squeeze=False)\n    info_by_slot={r[\"slot\"]:r for r in ci}\n    for ax,slot in zip(axes.flat,top):\n        ax.imshow(cx[slot,0],cmap=\"gray\"); ax.axis(\"off\")\n        r=info_by_slot[int(slot)]\n        ax.set_title(f\"{r['plane']} slice {r['slice_index']}\\nweight={w[slot]:.3f}\",fontsize=9)\n    savefig(\"10_attention_slices.png\")\n    print(\"高權重不等於已定位病灶；仍可能關注掃描特徵、偽影或其他線索。\")\n","id":"knee-035"},{"cell_type":"markdown","metadata":{},"source":"## 17｜補看預測的不穩定性\n\n我對同一個 OOF 模型啟用分類頭 dropout，重複預測 30 次，觀察分數變動。這是 MC dropout 的探索性指標，不是醫療可信區間。只在分類頭放 dropout 也只能反映部分模型變異；低標準差不能證明預測正確。\n","id":"knee-036"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    review_model.eval()\n    for module in review_model.modules():\n        if isinstance(module,nn.Dropout): module.train()\n    with torch.inference_mode():\n        draws=torch.stack([review_model(torch.as_tensor(X[case:case+1],device=DEVICE),\n                                       torch.as_tensor(M[case:case+1],device=DEVICE)).sigmoid()[0]\n                           for _ in range(30)]).cpu().numpy()\n    review_model.eval()\n    uncertainty=pd.DataFrame({\"target\":TARGETS,\"score_mean\":draws.mean(0),\"score_std\":draws.std(0),\"label\":Y[case]})\n    uncertainty.to_csv(OUT / \"case_mc_dropout.csv\",index=False); display(uncertainty)\n    plt.figure(figsize=(10,4)); plt.errorbar(range(12),draws.mean(0),yerr=draws.std(0),fmt=\"o\",capsize=3)\n    plt.xticks(range(12),TARGETS,rotation=45,ha=\"right\"); plt.ylim(0,1)\n    plt.ylabel(\"MC mean ± 1 standard deviation\"); savefig(\"11_mc_dropout.png\")\n","id":"knee-037"},{"cell_type":"markdown","metadata":{},"source":"## 18｜保存可獨立推論的模型\n\n我保存一份 CNN、五個 attention 分類頭、設定與解碼 wheel。正式推論不需要讀取報告、人工標籤或訓練影像。五個 fold 模型的輸出取平均，不再調整門檻，因為比賽要的是連續分數。\n","id":"knee-038"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    manifest=dict(config=CFG, targets=TARGETS, planes=PLANES, backbone_hash=BACKBONE_HASH,\n                  feature_signature=FEATURE_SIGNATURE, submission_pooling=\"attention\",\n                  trained_studies=len(gold), labeled_cells=int(np.isfinite(Y).sum()),\n                  versions=versions, checkpoints=[f\"attention_fold{f}.pt\" for f in range(CFG[\"folds\"])])\n    for name in manifest[\"checkpoints\"]: assert (OUT/name).exists()\n    (OUT/\"knee_manifest.json\").write_text(json.dumps(manifest,ensure_ascii=False,indent=2))\n    SOURCE=OUT\nelse:\n    assert manifest[\"backbone_hash\"]==BACKBONE_HASH, \"CNN 權重與 manifest 不一致。\"\n    assert manifest[\"feature_signature\"]==FEATURE_SIGNATURE\nprint(\"推論模型來源：\", SOURCE)\n","id":"knee-039"},{"cell_type":"markdown","metadata":{},"source":"## 19｜對測試資料產生提交檔\n\n我直接讀取當次 `test.csv` 的 study 清單。正式評分時 Kaggle 會把三筆示例換成完整隱藏資料，因此不能把那三個 UID 或 sample submission 的列數寫死。\n\n推論遇到解碼或序列問題會停止，不會靜默生成常數答案。每個 study 的耗時會保存，若預估接近 9 小時，我需先縮減切片或改善讀取，再重新訓練與保存一致設定。\n","id":"knee-040"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"test=pd.read_csv(ROOT/\"test.csv\",dtype={\"StudyInstanceUID\":str})\ntest_series=pd.read_csv(ROOT/\"test_series.csv\",dtype={\"StudyInstanceUID\":str,\"SeriesInstanceUID\":str})\nassert test.StudyInstanceUID.notna().all() and not test.StudyInstanceUID.duplicated().any()\nchosen=choose_series(test_series)\nseries_map={uid:g for uid,g in chosen.groupby(\"StudyInstanceUID\")}\nheads=[]\nfor name in manifest[\"checkpoints\"]:\n    model=KneeMIL(\"attention\").to(DEVICE)\n    ck=torch.load(SOURCE/name,map_location=\"cpu\",weights_only=True)\n    assert ck[\"pooling\"]==\"attention\" and ck[\"epoch\"]==CFG[\"epochs\"]\n    model.load_state_dict(ck[\"state\"]); model.eval(); heads.append(model)\n\npredictions=[]; timing=[]; inference_start=time.monotonic()\nfor uid in tqdm(test.StudyInstanceUID,desc=\"測試影像推論\"):\n    tick=time.monotonic()\n    feat,mask,_,_=get_features(uid,\"test\")\n    p=np.mean([predict_head(h,feat[None].astype(np.float32),mask[None].astype(bool))[0] for h in heads],axis=0)\n    predictions.append(p); timing.append(dict(StudyInstanceUID=uid,seconds=time.monotonic()-tick))\n    if len(timing)==3:\n        avg=np.mean([r[\"seconds\"] for r in timing])\n        print(f\"前三筆平均 {avg:.1f} 秒；1300 筆約 {avg*1300/3600:.2f} 小時（粗估，快取會低估）。\")\nsubmission=pd.DataFrame(np.stack(predictions),columns=TARGETS)\nsubmission.insert(0,\"StudyInstanceUID\",test.StudyInstanceUID.to_numpy())\nassert list(submission.columns)==[\"StudyInstanceUID\"]+TARGETS\nassert submission.StudyInstanceUID.tolist()==test.StudyInstanceUID.tolist()\nassert submission.shape==(len(test),13)\nassert np.isfinite(submission[TARGETS].to_numpy()).all()\nassert submission[TARGETS].ge(0).all().all() and submission[TARGETS].le(1).all().all()\nsubmission_path=Path(\"/kaggle/working/submission.csv\")\nsubmission.to_csv(submission_path,index=False)\npd.DataFrame(timing).to_csv(OUT/\"inference_timing.csv\",index=False)\ndisplay(submission.head())\nprint(\"已產生：\",submission_path,\"｜列數：\",len(submission),\"｜總推論秒數：\",round(time.monotonic()-inference_start,1))\n","id":"knee-041"},{"cell_type":"markdown","metadata":{},"source":"## 20｜把實驗結果整理成作品紀錄\n\n我只記錄實際跑出的結果，不預先寫「注意力一定比較好」。若注意力沒有勝過平均聚合，也是值得分析的結果：可能是標籤太少、固定 CNN 不適合 MRI，或目前取樣方式沒有保留足夠病灶資訊。\n\n輸出中的 OOF、folds 與案例包含競賽識別資訊；我保留作為私人實驗紀錄。公開作品時以程式、彙整數字與獲准使用的圖為主，影像、報告與衍生資料的分享須遵循競賽規範。\n","id":"knee-042"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"if MODE == \"train\":\n    a=summary.loc[\"attention\",\"macro_auc\"]; b=summary.loc[\"mean\",\"macro_auc\"]\n    valid=int(summary.loc[\"attention\",\"valid_targets\"])\n    text=f\"\"\"# 膝蓋 MRI 多標籤辨識：我的實驗紀錄\n\n## 動機\n因為自身膝蓋困擾，我希望透過 RSNA 比賽資料學習醫學影像處理，完成可重現的研究基準。\n\n## 方法\n三切面各最多 {CFG['slices_per_plane']} 張影像、固定 ImageNet ResNet-18、平均／注意力聚合、12 類 masked BCE。\n使用 {len(gold)} 筆有人工標籤的檢查與 {CFG['folds']} folds，固定 {CFG['epochs']} epochs。\n有 PatientID 時以其分組，否則以 study 分組；未完成院所外部驗證。\n報告僅用於探索，未做多模態訓練或弱標籤擴增。\n\n## 實測 OOF 結果\n- mean 平均 AUC：{b:.4f}\n- attention 平均 AUC：{a:.4f}\n- 差異：{a-b:+.4f}\n- 可評估類別：{valid}/12\n\n這是本地 OOF 結果，並非 Kaggle leaderboard 分數。樣本少且未有獨立外部測試，不能宣稱臨床有效。\n正式提交結果：尚待 Kaggle 評分。\n\n## 我需要進一步檢查的問題\n1. 注意力有沒有穩定勝過平均聚合，而不只是單次切分差異？\n2. 模型是否依賴儀器、序列或影像邊緣線索？\n3. 哪些類別陽性樣本最少，是否主導分數變動？\n4. 多語言報告的否定、程度與不確定語句，如何轉成經核對的弱標籤？\n\n## 技術貢獻與界線\n本專案整合 DICOM 幾何排序、解碼、跨切面取樣、固定 CNN 特徵、MIL 聚合、分組驗證、可續跑 checkpoint 與離線推論。\nResNet、注意力 MIL 與 MC dropout 均為既有方法；本階段貢獻是應用與實驗整合，不宣稱發明新的演算法。\n\"\"\"\n    (OUT/\"project_summary.md\").write_text(text,encoding=\"utf-8\")\n    display(Markdown(text))\n","id":"knee-043"},{"cell_type":"markdown","metadata":{},"source":"## 21｜下一步我想研究的方向\n\n第一版完成的是端到端基準管線，並不等於已把整個比賽資料充分利用。我會依結果選擇下一步，而不是一次增加所有方法。\n\n1. **弱監督標籤**：用在 Kaggle 本地執行的多語言模型，從報告抽取各類的「存在／不存在／不確定／未提及」與證據片段。未提及不當成陰性，並檢查官方嚴重度定義；保留人工標籤作獨立驗證。正式測試仍只輸入影像。\n2. **領域適應**：比較 ImageNet、MRI 自監督特徵與微調 CNN。若做預訓練，需排除外部驗證病人，並記錄外部資料來源。\n3. **序列與切片策略**：固定總切片數比較一個／三個切面，再研究各疾病不同的注意力。這樣比較才不會把切片數增加誤認成架構改善。\n4. **更可靠的評估**：增加人工核對標籤、使用病人與院所分組，並在新來源資料做外部驗證。小樣本 OOF 不能替代這些工作。\n\n推甄時，我可以先用「**結合多切面 MRI 特徵與注意力聚合之膝關節異常辨識實驗**」作為題名，說明自己的動機、實作選擇、對照實驗與尚未解決的問題。等跑出結果並理解程式後，再把它寫成已完成的作品成果。\n","id":"knee-044"},{"cell_type":"markdown","metadata":{},"source":"## 22｜輸出檔案與我會留下的證據\n\n| 檔案 | 用途 |\n|---|---|\n| `submission.csv`（working 根目錄） | Kaggle 提交 |\n| `knee_project/backbone.pt` | 固定 CNN 權重 |\n| `attention_fold*.pt` / `mean_fold*.pt` | 分類頭、optimizer、續跑狀態 |\n| `knee_manifest.json` / `environment.json` | 設定與版本 |\n| `wheels/` | 提交時離線安裝 DICOM 解碼依賴 |\n| `features/` | 每個 study 的小型特徵快取 |\n| `folds.csv` / `oof_*.csv` | 分組切分與真正 OOF 預測 |\n| `metrics_by_target.csv` / `experiment_summary.csv` | 逐類與整體結果 |\n| `history.csv` / `bootstrap.csv` | 訓練歷程與抽樣分析 |\n| `figures/` / `error_cases.csv` | 圖表與錯誤案例 |\n| `project_summary.md` | 根據實測結果產生的作品紀錄 |\n\n我會保留訓練 Notebook 的私人 Output，提交版本掛載整個輸出即可，不需要手動逐一下載大型 MRI。\n","id":"knee-045"},{"cell_type":"code","metadata":{},"execution_count":null,"outputs":[],"source":"files=[dict(file=str(p.relative_to(OUT)),MB=round(p.stat().st_size/1024**2,3))\n       for p in OUT.rglob(\"*\") if p.is_file() and p.parent.name!=\"features\"]\ndisplay(pd.DataFrame(files).sort_values(\"file\").reset_index(drop=True))\nprint(\"目前輸出總大小（GB）：\",round(sum(p.stat().st_size for p in OUT.rglob(\"*\") if p.is_file())/1024**3,3))\nprint(\"完成後我還需要確認 Kaggle 版本保存成功；正式提交需要另外的 infer 版本。\")\n","id":"knee-046"},{"cell_type":"markdown","metadata":{},"source":"## 23｜資料來源與方法參考\n\n比賽頁已於 2026-09-05 核對；Notebook 會再依實際 CSV 檢查欄位與標籤數。\n\n- [Kaggle：官方任務、評分、期限與執行限制](https://www.kaggle.com/competitions/rsna-knee-abnormality-detection/overview)\n- [Kaggle：資料欄位、DICOM 格式與測試資料說明](https://www.kaggle.com/competitions/rsna-knee-abnormality-detection/data)\n- [RSNA：Knee Abnormality Detection AI Challenge](https://www.rsna.org/artificial-intelligence/ai-image-challenge/knee-mri-ai-challenge)\n- [PyTorch：ResNet-18 與預訓練權重](https://docs.pytorch.org/vision/stable/models/generated/torchvision.models.resnet18.html)\n- [pydicom：壓縮影像解碼](https://pydicom.github.io/pydicom/stable/guides/user/image_data_handlers.html)\n- [Ilse et al., *Attention-based Deep Multiple Instance Learning*, ICML 2018](https://proceedings.mlr.press/v80/ilse18a)。注意力聚合的概念來源；這份程式是簡化實作。\n- [Gal & Ghahramani, *Dropout as a Bayesian Approximation*, ICML 2016](https://proceedings.mlr.press/v48/gal16.pdf)。MC dropout 的概念來源。\n\n準備時也透過 GitHub 核對公開參賽者的資料整理；最終欄位與提交規格以上述官方頁為準。本筆記未複製其訓練程式、標籤或權重：\n[competition reference](https://github.com/Beiciccc/rsna-knee-abnormality-detection/blob/main/docs/competition.md)。\n","id":"knee-047"}],"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"name":"python","version":"3.11.0"},"kaggle":{"title":"膝蓋 MRI 異常辨識｜正式提交版"}},"nbformat":4,"nbformat_minor":5}