{"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.11.13"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":99552,"databundleVersionId":13851420,"sourceType":"competition"},{"sourceId":12746899,"sourceType":"datasetVersion","datasetId":8057988}],"dockerImageVersionId":31090,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Configs and Import Libraries","metadata":{}},{"cell_type":"code","source":"# ==== Imports ====\nfrom tqdm.auto import tqdm\nimport os\nimport shutil\nimport numpy as np\nimport pandas as pd\nimport pydicom\nfrom sklearn.model_selection import StratifiedKFold\nfrom sklearn.metrics import roc_auc_score\n\nimport torch\nimport torch.nn as nn\nfrom torch.utils.data import Dataset, DataLoader\nimport timm\nimport albumentations as A\nfrom albumentations.pytorch import ToTensorV2\nimport cv2\nimport gc\nimport ast\nimport sys\n\nimport matplotlib.pyplot as plt\nfrom matplotlib import patches\n\nfrom torch.utils.data import Dataset\nimport torch\n\n# ==== Config ====\nUSE_CV = True\nNUM_FOLDS = 5\nEPOCHS = 5\nBATCH_SIZE = 8\nLR = 1e-4\nIMG_SIZE = 224 #224\nSERIES_ROOT = \"/kaggle/input/rsna-intracranial-aneurysm-detection/series\"\nTRAIN_CSV = \"/kaggle/input/rsna-intracranial-aneurysm-detection/train.csv\"\nCACHE_DIR = \"/kaggle/input/rsna-2025-ia-cta-224-tensors\"\nwindowing = False\n\nARCHITECTURES = ['tf_efficientnet_b0'] #,'resnet18'\n\nLABEL_COLS = [\n    'Left Infraclinoid Internal Carotid Artery', 'Right Infraclinoid Internal Carotid Artery',\n    'Left Supraclinoid Internal Carotid Artery', 'Right Supraclinoid Internal Carotid Artery',\n    'Left Middle Cerebral Artery', 'Right Middle Cerebral Artery',\n    'Anterior Communicating Artery', 'Left Anterior Cerebral Artery',\n    'Right Anterior Cerebral Artery', 'Left Posterior Communicating Artery',\n    'Right Posterior Communicating Artery', 'Basilar Tip',\n    'Other Posterior Circulation', 'Aneurysm Present'\n]\nID_COL = \"SeriesInstanceUID\"\n\n","metadata":{"execution":{"iopub.status.busy":"2025-10-12T18:37:26.948333Z","iopub.execute_input":"2025-10-12T18:37:26.94861Z","iopub.status.idle":"2025-10-12T18:37:40.875363Z","shell.execute_reply.started":"2025-10-12T18:37:26.948588Z","shell.execute_reply":"2025-10-12T18:37:40.874587Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Use GPU or CPU and related Settings","metadata":{}},{"cell_type":"code","source":"# ==== Detect GPUs ====\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nnum_gpus = torch.cuda.device_count()\nUSE_MULTIGPU = num_gpus > 1\nprint(f\"Detected {num_gpus} GPU(s)\")\n\n# ==== Set up CPU/GPU ====\nimport os, torch\n\nNUM_CPUS = os.cpu_count() or 2\n# Use most cores for compute; keep a couple free for I/O/OS\ntorch.set_num_threads(max(1, NUM_CPUS - 2))\ntorch.set_num_interop_threads(2)\n\n\nNUM_WORKERS = max(4, min(8, NUM_CPUS // 2))  # tune as needed\n\nos.environ[\"OMP_NUM_THREADS\"] = str(max(1, NUM_CPUS - 2))\nos.environ[\"MKL_NUM_THREADS\"] = str(max(1, NUM_CPUS - 2))\n\nprint(\"CPUs:\", NUM_CPUS)\nprint(\"Torch threads:\", torch.get_num_threads(), \"interop:\", torch.get_num_interop_threads())\n\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-12T18:37:40.876543Z","iopub.execute_input":"2025-10-12T18:37:40.876818Z","iopub.status.idle":"2025-10-12T18:37:40.964399Z","shell.execute_reply.started":"2025-10-12T18:37:40.876789Z","shell.execute_reply":"2025-10-12T18:37:40.963576Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Define some useful functions","metadata":{}},{"cell_type":"code","source":"# ==== Utils ====\ndef get_model(arch):\n    model = timm.create_model(arch, in_chans=IN_CHANS, num_classes=len(LABEL_COLS), pretrained=False)\n    if USE_MULTIGPU:\n        model = nn.DataParallel(model)\n    return model.to(device)\n\ndef sort_dicom_slices(filepaths):\n    dicoms = [pydicom.dcmread(fp, force=True) for fp in filepaths]\n    try:\n        dicoms.sort(key=lambda d: float(d.ImagePositionPatient[2]))\n    except Exception:\n        dicoms.sort(key=lambda d: int(getattr(d, \"InstanceNumber\", 0)))\n    return dicoms\n\ndef to_hu(dcm):\n    img = dcm.pixel_array.astype(np.float32)\n    slope = float(getattr(dcm, \"RescaleSlope\", 1.0))\n    intercept = float(getattr(dcm, \"RescaleIntercept\", 0.0))\n    return img * slope + intercept\n\ndef window_img(hu, wl, ww):\n    lo, hi = wl - ww/2.0, wl + ww/2.0\n    hu = np.clip(hu, lo, hi)\n    return ((hu - lo) / (hi - lo)).astype(np.float32)\n\n# Helper function to parse the coordinates from the format {x: 12.213, y: 342.132}\ndef parse_coordinates(coord_str):\n    coord_dict = ast.literal_eval(coord_str)\n    return np.array([coord_dict['x'], coord_dict['y']])\n\n# Load the train_localizers.csv\nlocalizer_df = pd.read_csv('/kaggle/input/rsna-intracranial-aneurysm-detection/train_localizers.csv')\n\n# Apply the parse function to the coordinates column\nlocalizer_df['coordinates'] = localizer_df['coordinates'].apply(parse_coordinates)\n\n\ndef valid_coords(coords):\n    return np.all(np.isfinite(coords)) and not (coords[0] == 0.0 and coords[1] == 0.0)\n\ndef coords_to_px(coords, img_size=IMG_SIZE):\n    x, y = float(coords[0]), float(coords[1])\n    if 0.0 <= x <= 1.0 and 0.0 <= y <= 1.0:\n        x *= img_size; y *= img_size\n    return int(round(x)), int(round(y))\n\ndef make_bbox_px(x, y, img_size=IMG_SIZE, box_frac=0.15, min_px=24, max_px=None):\n    r = max(min_px/2, box_frac * img_size / 2.0)\n    if max_px is not None:\n        r = min(r, max_px/2)\n    x1 = int(np.clip(x - r, 0, img_size - 1))\n    y1 = int(np.clip(y - r, 0, img_size - 1))\n    x2 = int(np.clip(x + r, 0, img_size - 1))\n    y2 = int(np.clip(y + r, 0, img_size - 1))\n    # enforce minimum width/height of 2 px\n    if x2 <= x1: x2 = min(img_size - 1, x1 + 1)\n    if y2 <= y1: y2 = min(img_size - 1, y1 + 1)\n    return x1, y1, x2, y2\n\ndef crop_and_resize_chw(img_chw, x1, y1, x2, y2, out_size=IMG_SIZE):\n    # img_chw: [C,H,W] float16/float32\n    img_hwc = np.transpose(np.asarray(img_chw), (1, 2, 0))\n    crop = img_hwc[y1:y2, x1:x2]\n    if crop.size == 0 or crop.shape[0] < 2 or crop.shape[1] < 2:\n        # fallback to full image (cast to float32 to avoid later issues)\n        full = img_hwc.astype(np.float32, copy=False)\n        full = cv2.resize(full, (out_size, out_size), interpolation=cv2.INTER_AREA)\n        return np.transpose(full, (2, 0, 1))\n\n    crop = crop.astype(np.float32, copy=False)\n    crop = np.ascontiguousarray(crop)\n    crop = cv2.resize(crop, (out_size, out_size), interpolation=cv2.INTER_AREA)\n    return np.transpose(crop, (2, 0, 1))\n\ndef crop_and_resize_hwc(img_hwc: np.ndarray, x1: int, y1: int, x2: int, y2: int, out_size: int) -> np.ndarray:\n    crop = img_hwc[y1:y2, x1:x2]\n    if crop.size == 0 or crop.shape[0] < 2 or crop.shape[1] < 2:\n        crop = img_hwc  # fallback to full img\n    crop = crop.astype(np.float32, copy=False)\n    crop = np.ascontiguousarray(crop)\n    crop = cv2.resize(crop, (out_size, out_size), interpolation=cv2.INTER_AREA)\n    return crop  # HWC float32\n\n\ndef make_loader(ds, shuffle, num_workers):\n    kwargs = dict(batch_size=BATCH_SIZE, shuffle=shuffle, pin_memory=True)\n    if num_workers > 0:\n        kwargs.update(num_workers=num_workers, persistent_workers=False, prefetch_factor=2)\n    return DataLoader(ds, **kwargs)\n\n# In notebook, prefer workers=0 for clean shutdown\nfrom IPython import get_ipython\nIN_NOTEBOOK = get_ipython() is not None\nSAFE_WORKERS = 0 if IN_NOTEBOOK else NUM_WORKERS\n\n\n\ndef compute_weighted_auc(y_true, y_prob, class_names):\n    y_true = np.atleast_2d(y_true)\n    y_prob = np.atleast_2d(y_prob)\n    aucs, skipped = {}, []\n    for i, name in enumerate(class_names):\n        yi = y_true[:, i]\n        if len(np.unique(yi)) < 2:\n            skipped.append(name)\n            continue\n        aucs[name] = roc_auc_score(yi, y_prob[:, i])\n    ap_name = \"Aneurysm Present\"\n    ap_auc = aucs.get(ap_name, np.nan)\n    others = [v for k, v in aucs.items() if k != ap_name]\n    others_mean = np.mean(others) if others else np.nan\n    weighted_auc = 0.5 * (ap_auc + others_mean)\n    return weighted_auc, ap_auc, others_mean, aucs, skipped\n\n\n# 2.5D window config\nOFFSETS = (-2, -1, 0, 1, 2)  # window length = 5 (set to (-1,0,1) for 3)\nWINDOW_LEN = len(OFFSETS)\n\ndef load_cached_img(img_path: str) -> np.ndarray:\n    if img_path.endswith(\".npz\"):\n        return np.load(img_path)[\"arr_0\"]          # [C_all, H, W]\n    return np.load(img_path, mmap_mode=\"r\")        # [C_all, H, W]\n\ndef take_window(img_chw: np.ndarray, center_idx: int, offsets=OFFSETS) -> np.ndarray:\n    C = img_chw.shape[0]\n    idxs = [min(max(0, center_idx + o), C - 1) for o in offsets]\n    return img_chw[idxs, :, :]                     # [WINDOW_LEN, H, W]\n\n\n@torch.no_grad()\ndef infer_series_all_slices(model: nn.Module, sid: str, cache_dir=CACHE_DIR, img_size=IMG_SIZE,\n                            offsets=OFFSETS, batch_size=16, aggregate=\"max\", use_roi=True):\n    img_path_npz = os.path.join(cache_dir, f\"{sid}_img.npz\")\n    img_path_npy = os.path.join(cache_dir, f\"{sid}_img.npy\")\n    img_path = img_path_npz if os.path.exists(img_path_npz) else img_path_npy\n\n    coord_path = os.path.join(cache_dir, f\"{sid}_coords.npy\")\n\n    img_full = load_cached_img(img_path)                  # [C_all, H, W]\n    coords = np.load(coord_path).astype(np.float32) if os.path.exists(coord_path) else np.array([0., 0.], np.float32)\n    C_all = img_full.shape[0]\n\n    def prep_window(center_idx: int) -> tuple[np.ndarray, np.ndarray]:\n        win = take_window(img_full, center_idx, offsets)   # [WINDOW_LEN, H, W]\n        if use_roi and valid_coords(coords):\n            cx, cy = coords_to_px(coords, img_size)\n            x1, y1, x2, y2 = make_bbox_px(cx, cy, img_size, 0.15, 24)\n            win_roi = crop_and_resize_chw(win, x1, y1, x2, y2, img_size)  # [WINDOW_LEN,H,W]\n        else:\n            win_hwc = np.transpose(win, (1, 2, 0)).astype(np.float32, copy=False)\n            win_hwc = np.ascontiguousarray(win_hwc)\n            resized = cv2.resize(win_hwc, (img_size, img_size), interpolation=cv2.INTER_AREA)\n            win_roi = np.transpose(resized, (2, 0, 1))\n        win_full = win_roi  # if you want two branches identical spatially; else keep uncropped for full\n        return win_full, win_roi\n\n    xs_full, xs_roi = [], []\n    for c in range(C_all):\n        w_full, w_roi = prep_window(c)\n        xs_full.append(w_full)\n        xs_roi.append(w_roi)\n\n    xs_full = torch.from_numpy(np.stack(xs_full).astype(np.float32)).to(device)  # [N, C, H, W]\n    xs_roi  = torch.from_numpy(np.stack(xs_roi).astype(np.float32)).to(device)\n\n    logits_all = []\n    for i in range(0, xs_full.shape[0], batch_size):\n        xb_full = xs_full[i:i+batch_size]\n        xb_roi  = xs_roi[i:i+batch_size]\n        coords_b = torch.from_numpy(np.repeat(coords[None, :], xb_full.shape[0], axis=0)).to(device)\n        logits = model(xb_full, xb_roi, coords_b)  # [B, num_classes]\n        logits_all.append(logits)\n    logits_all = torch.cat(logits_all, dim=0)\n    probs_all = torch.sigmoid(logits_all).cpu().numpy()     # [N_windows, num_classes]\n\n    if aggregate == \"max\":\n        series_prob = probs_all.max(axis=0)\n    elif aggregate == \"mean\":\n        series_prob = probs_all.mean(axis=0)\n    elif aggregate == \"topk_mean\":\n        k = max(1, xs_full.shape[0] // 5)\n        series_prob = np.sort(probs_all, axis=0)[-k:].mean(axis=0)\n    else:\n        series_prob = probs_all.mean(axis=0)\n    return series_prob\n","metadata":{"execution":{"iopub.status.busy":"2025-10-12T18:37:40.965124Z","iopub.execute_input":"2025-10-12T18:37:40.965351Z","iopub.status.idle":"2025-10-12T18:37:41.055451Z","shell.execute_reply.started":"2025-10-12T18:37:40.965334Z","shell.execute_reply":"2025-10-12T18:37:41.054918Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Preview Slices and Region of Interest","metadata":{}},{"cell_type":"code","source":"# ==== Preview full image and ROI from cache ====\n\n\ndef show_cached_full_and_roi(\n    sid: str,\n    cache_dir: str = CACHE_DIR,\n    img_size: int = IMG_SIZE,\n    roi_box_frac: float = 0.15,\n    roi_min_px: int = 24,\n):\n    # Load cached CHW and coords (support .npy and .npz)\n    img_path_npz = os.path.join(cache_dir, f\"{sid}_img.npz\")\n    img_path_npy = os.path.join(cache_dir, f\"{sid}_img.npy\")\n    img_path = img_path_npz if os.path.exists(img_path_npz) else img_path_npy\n    img_chw = load_cached_img(img_path)  # [C,H,W]\n    coords_path = os.path.join(cache_dir, f\"{sid}_coords.npy\")\n    coords = np.load(coords_path).astype(np.float32) if os.path.exists(coords_path) else np.array([0.0, 0.0], dtype=np.float32)\n\n    # CHW -> HWC\n    img_hwc = np.transpose(np.asarray(img_chw), (1, 2, 0)).astype(np.float32, copy=False)\n\n    # Visualization image (grayscale): mean over channels\n    vis_full = img_hwc.mean(axis=2)\n\n    # Compute ROI box using actual image size\n    H, W = img_hwc.shape[:2]\n    side = min(H, W)\n    if valid_coords(coords):\n        cx, cy = coords_to_px(coords, side)\n        x1, y1, x2, y2 = make_bbox_px(cx, cy, side, roi_box_frac, roi_min_px)\n        vis_roi_hwc = crop_and_resize_hwc(img_hwc, x1, y1, x2, y2, img_size)\n        vis_roi = vis_roi_hwc.mean(axis=2)\n    else:\n        # Fallback: show full image as ROI if coords invalid\n        x1 = y1 = 0\n        x2 = W - 1\n        y2 = H - 1\n        vis_roi = vis_full\n\n    # Plot\n    fig, axs = plt.subplots(1, 2, figsize=(10, 5))\n    axs[0].imshow(vis_full, cmap=\"gray\")\n    axs[0].set_title(f\"Full ({img_chw.shape[0]}-slice) 2.5D\")\n    rect = patches.Rectangle(\n        (x1, y1), x2 - x1, y2 - y1, linewidth=2, edgecolor=\"lime\", facecolor=\"none\"\n    )\n    axs[0].add_patch(rect)\n    axs[0].axis(\"off\")\n\n    axs[1].imshow(vis_roi, cmap=\"gray\")\n    axs[1].set_title(\"ROI (resized)\")\n    axs[1].axis(\"off\")\n\n    plt.tight_layout()\n    plt.show()\n\n# Example usage:\n# Pick a SeriesInstanceUID present in your cache\ndf = pd.read_csv(TRAIN_CSV)\ndf = df[df['Modality'] == 'CTA']\n\nexample_sid = str(df.iloc[0][ID_COL])\nshow_cached_full_and_roi(example_sid)","metadata":{"execution":{"iopub.status.busy":"2025-10-12T18:37:41.056962Z","iopub.execute_input":"2025-10-12T18:37:41.05719Z","iopub.status.idle":"2025-10-12T18:37:42.28374Z","shell.execute_reply.started":"2025-10-12T18:37:41.057172Z","shell.execute_reply":"2025-10-12T18:37:42.283013Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Random Showing of slices from a series","metadata":{}},{"cell_type":"code","source":"\nsid = os.path.basename(np.random.choice([p for p in os.listdir(CACHE_DIR) if p.endswith('_img.npy') or p.endswith('_img.npz')])).replace('_img.npy','').replace('_img.npz','')\nprint(sid)\nshow_cached_full_and_roi(sid)","metadata":{"execution":{"iopub.status.busy":"2025-10-12T18:37:42.284461Z","iopub.execute_input":"2025-10-12T18:37:42.284687Z","iopub.status.idle":"2025-10-12T18:37:43.07905Z","shell.execute_reply.started":"2025-10-12T18:37:42.284668Z","shell.execute_reply":"2025-10-12T18:37:43.07828Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Example Series with Slices that do not have ROI","metadata":{}},{"cell_type":"code","source":"#no ROI\nshow_cached_full_and_roi(\"1.2.826.0.1.3680043.8.498.99348616762869189162377000017558094817\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-12T18:37:43.079794Z","iopub.execute_input":"2025-10-12T18:37:43.080028Z","iopub.status.idle":"2025-10-12T18:37:44.70117Z","shell.execute_reply.started":"2025-10-12T18:37:43.080009Z","shell.execute_reply":"2025-10-12T18:37:44.700391Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def debug_roi_bbox(sid: str, cache_dir: str = CACHE_DIR):\n    import numpy as np, os\n    img_path_npz = os.path.join(cache_dir, f\"{sid}_img.npz\")\n    img_path_npy = os.path.join(cache_dir, f\"{sid}_img.npy\")\n    img_path = img_path_npz if os.path.exists(img_path_npz) else img_path_npy\n    img_chw = load_cached_img(img_path)\n    coords_path = os.path.join(cache_dir, f\"{sid}_coords.npy\")\n    coords = np.load(coords_path).astype(np.float32) if os.path.exists(coords_path) else np.array([0.0, 0.0], dtype=np.float32)\n    print(\"img shape (C,H,W):\", img_chw.shape, \"coords:\", coords)\n\n    H = W = img_chw.shape[1]\n    cx, cy = coords_to_px(coords, W)\n    x1, y1, x2, y2 = make_bbox_px(cx, cy, W, box_frac=0.15, min_px=24)\n    print(\"px center:\", (cx, cy), \"bbox:\", (x1, y1, x2, y2))\n\n    if (x2 - x1) < 2 or (y2 - y1) < 2:\n        print(\"ROI degenerated → fallback to full image\")\n    if not valid_coords(coords):\n        print(\"coords invalid per valid_coords → ROI fallback to full image\")","metadata":{"execution":{"iopub.status.busy":"2025-10-12T18:37:44.702142Z","iopub.execute_input":"2025-10-12T18:37:44.702664Z","iopub.status.idle":"2025-10-12T18:37:44.708817Z","shell.execute_reply.started":"2025-10-12T18:37:44.702637Z","shell.execute_reply":"2025-10-12T18:37:44.708223Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Model Definition","metadata":{}},{"cell_type":"code","source":"\n# Minimal cached dataset (uses mmap for low-RAM I/O)\n\n# Example transforms (no ToTensorV2; we convert to tensor ourselves)\ntrain_tf = A.Compose([\n    A.HorizontalFlip(p=0.5),\n    A.Affine(\n        scale=(0.95,1.05),\n        translate_percent={'x':(-0.05,0.05),'y':(-0.05,0.05)},\n        rotate=(-10,10),\n        p=0.5\n    ),\n    A.GaussianNoise(var_limit=(1.0, 5.0), mean=0, p=0.2) if hasattr(A, \"GaussianNoise\") else A.GaussNoise(p=0.2),\n    A.MotionBlur(blur_limit=3, p=0.1),\n])\nval_tf = None  # keep cached normalization; or add light Normalize if desired\n\n# ===== 2) Hybrid (full + ROI) Cached Dataset with Albumentations =====\nclass HybridCachedNPYDatasetAlb(Dataset):\n    def __init__(self, df: pd.DataFrame, id_col: str, img_size: int,\n                 roi_box_frac: float = 0.2, roi_min_px: int = 10,\n                 transform: A.Compose | None = None):\n        # Filter to SIDs that have cached image files\n        all_files = []\n        try:\n            all_files = os.listdir(CACHE_DIR)\n        except Exception:\n            all_files = []\n        avail = ({f.replace('_img.npy', '') for f in all_files if f.endswith('_img.npy')} |\n                 {f.replace('_img.npz', '') for f in all_files if f.endswith('_img.npz')})\n        df = df[df[id_col].astype(str).isin(avail)].reset_index(drop=True)\n\n        self.df = df\n        self.id_col = id_col\n        self.img_size = img_size\n        self.roi_box_frac = roi_box_frac\n        self.roi_min_px = roi_min_px\n        # Use additional_targets to apply the same augmentation to both full and ROI\n        if transform is None:\n            self.transform = A.Compose(\n                [\n                    A.HorizontalFlip(p=0.5),\n                    A.Affine(\n                        scale=(0.95,1.05),\n                        translate_percent={'x':(-0.05,0.05),'y':(-0.05,0.05)},\n                        rotate=(-10,10),\n                        p=0.5,\n                    ),\n                ],\n                additional_targets={\"image2\": \"image\"},\n            )\n        else:\n            # Ensure additional_targets exists\n            self.transform = A.Compose(transform.transforms, additional_targets={\"image2\": \"image\"})\n\n    def __len__(self):\n        return len(self.df)\n\n    def __getitem__(self, idx: int):\n        row = self.df.iloc[idx]\n        sid = str(row[self.id_col])\n\n        img_path_npz = os.path.join(CACHE_DIR, f\"{sid}_img.npz\")\n        img_path_npy = os.path.join(CACHE_DIR, f\"{sid}_img.npy\")\n        img_path = img_path_npz if os.path.exists(img_path_npz) else img_path_npy\n        coord_path = os.path.join(CACHE_DIR, f\"{sid}_coords.npy\")\n\n        img_chw = load_cached_img(img_path)              # [C,H,W]\n        # Robust to missing coords file\n        try:\n            coords = np.load(coord_path).astype(np.float32)     # [2]\n        except FileNotFoundError:\n            coords = np.array([0.0, 0.0], dtype=np.float32)\n\n        # Take a fixed window of slices to produce consistent channels (WINDOW_LEN)\n        C_all = img_chw.shape[0]\n        center = C_all // 2\n        img_win = take_window(img_chw, center, OFFSETS)  # [WINDOW_LEN, H, W]\n\n        # Full image HWC from window\n        img_full_hwc = np.transpose(np.asarray(img_win), (1, 2, 0)).astype(np.float32, copy=False)\n\n        # ROI image HWC from window\n        if valid_coords(coords):\n            cx, cy = coords_to_px(coords, self.img_size)\n            x1, y1, x2, y2 = make_bbox_px(cx, cy, self.img_size, self.roi_box_frac, self.roi_min_px)\n            img_roi_hwc = crop_and_resize_hwc(img_full_hwc, x1, y1, x2, y2, self.img_size)\n        else:\n            img_roi_hwc = img_full_hwc\n\n        # Ensure same H,W for both streams before augmentation\n        img_full_hwc = img_full_hwc.astype(np.float32, copy=False)\n        img_roi_hwc  = img_roi_hwc.astype(np.float32, copy=False)\n        img_full_hwc = np.ascontiguousarray(img_full_hwc)\n        img_roi_hwc  = np.ascontiguousarray(img_roi_hwc)\n        img_full_hwc = cv2.resize(img_full_hwc, (self.img_size, self.img_size), interpolation=cv2.INTER_AREA)\n        img_roi_hwc  = cv2.resize(img_roi_hwc,  (self.img_size, self.img_size), interpolation=cv2.INTER_AREA)\n\n        # Augment both consistently\n        aug = self.transform(image=img_full_hwc, image2=img_roi_hwc)\n        img_full_hwc = aug[\"image\"]\n        img_roi_hwc = aug[\"image2\"]\n\n        # Back to CHW and tensors\n        x_full = torch.from_numpy(np.transpose(img_full_hwc, (2, 0, 1)).copy()).float()\n        x_roi  = torch.from_numpy(np.transpose(img_roi_hwc,  (2, 0, 1)).copy()).float()\n\n        vals = pd.to_numeric(row[LABEL_COLS], errors=\"coerce\").values.astype(np.float32)\n        vals = np.nan_to_num(vals, nan=0.0)\n        y = torch.from_numpy(vals)\n\n        c = torch.from_numpy(coords)\n        return x_full, x_roi, y, c","metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","execution":{"iopub.status.busy":"2025-10-12T18:37:44.709547Z","iopub.execute_input":"2025-10-12T18:37:44.710296Z","iopub.status.idle":"2025-10-12T18:37:44.728989Z","shell.execute_reply.started":"2025-10-12T18:37:44.710268Z","shell.execute_reply":"2025-10-12T18:37:44.728246Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\n\nclass UNetBackbone(nn.Module):\n    \"\"\"U-Net编码器，支持通过参数调整网络深度\"\"\"\n    def __init__(self, in_chans=5, base_chans=64, depth=4):\n        super().__init__()\n        layers = []\n        current_chans = in_chans\n        \n        # 根据深度构建不同复杂度的U-Net\n        for i in range(depth):\n            # 通道数翻倍\n            next_chans = base_chans * (2 ** i) if i > 0 else base_chans\n            layers.extend([\n                nn.Conv2d(current_chans, next_chans, kernel_size=3, padding=1),\n                nn.BatchNorm2d(next_chans),\n                nn.ReLU(inplace=True),\n                nn.Conv2d(next_chans, next_chans, kernel_size=3, padding=1),\n                nn.BatchNorm2d(next_chans),\n                nn.ReLU(inplace=True),\n                nn.MaxPool2d(2)  # 下采样\n            ])\n            current_chans = next_chans\n        \n        # 最终层固定为1280通道，匹配原始模型\n        layers.extend([\n            nn.Conv2d(current_chans, 1280, kernel_size=3, padding=1),\n            nn.BatchNorm2d(1280),\n            nn.ReLU(inplace=True)\n        ])\n        \n        self.encoder = nn.Sequential(*layers)\n        self.feature_info = [{'num_chs': 1280}]\n\n    def forward(self, x):\n        return [self.encoder(x)]\n\n\nclass HybridAneurysmModel(nn.Module):\n    # 保留base_model_name参数以兼容训练代码，但实际使用U-Net\n    def __init__(self, base_model_name=\"unet_backbone\", num_classes=14):\n        super().__init__()\n\n        # 可以根据base_model_name调整U-Net的深度\n        depth_map = {\n            \"tf_efficientnet_attn_b5\": 4,  # 对应原始模型的复杂度\n            \"unet_small\": 3,\n            \"unet_large\": 5\n        }\n        depth = depth_map.get(base_model_name, 4)  # 默认深度\n        \n        # 1. U-Net编码器（兼容原始in_chans=5）\n        self.backbone = UNetBackbone(\n            in_chans=5,\n            depth=depth\n        )\n\n        # 2. 最终层通道数（与原始逻辑一致）\n        enc_dim = self.backbone.feature_info[-1]['num_chs']\n\n        # 3. 坐标小网络（与原始完全一致）\n        self.coord_fc = nn.Sequential(\n            nn.Linear(2, 32), nn.ReLU(), nn.Linear(32, 64)\n        )\n\n        # 4. 分类头（与原始完全一致）\n        self.classifier = nn.Sequential(\n            nn.Dropout(0.3),\n            nn.Linear(enc_dim * 2 + 64, num_classes)\n        )\n\n    def forward(self, x_full: torch.Tensor, x_roi: torch.Tensor, coords: torch.Tensor):\n        # 与原始模型前向传播完全一致\n        feat_full = self.backbone(x_full)[-1]\n        feat_roi  = self.backbone(x_roi)[-1]\n\n        att_full = torch.softmax(feat_full.view(*feat_full.shape[:2], -1), dim=-1) \\\n                          .view_as(feat_full)\n        att_roi  = torch.softmax(feat_roi.view(*feat_roi.shape[:2], -1), dim=-1) \\\n                          .view_as(feat_roi)\n\n        g_full = (feat_full * att_full).sum(dim=(2, 3))\n        g_roi  = (feat_roi  * att_roi).sum(dim=(2, 3))\n\n        f_coord = self.coord_fc(coords.float())\n\n        fused = torch.cat([g_full, g_roi, f_coord], dim=1)\n        return self.classifier(fused)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-12T18:37:44.729769Z","iopub.execute_input":"2025-10-12T18:37:44.73053Z","iopub.status.idle":"2025-10-12T18:37:44.74849Z","shell.execute_reply.started":"2025-10-12T18:37:44.730511Z","shell.execute_reply":"2025-10-12T18:37:44.747831Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Model Training Configurations","metadata":{}},{"cell_type":"code","source":"\n\ndf = pd.read_csv(TRAIN_CSV)\ndf = df[df['Modality'] == 'CTA']\n\n# Keep only SIDs present in cache\ntry:\n    files = os.listdir(CACHE_DIR)\n    cache_ids_npy = {f.replace('_img.npy','') for f in files if f.endswith('_img.npy')}\n    cache_ids_npz = {f.replace('_img.npz','') for f in files if f.endswith('_img.npz')}\n    cache_ids = cache_ids_npy | cache_ids_npz\n    df = df[df[ID_COL].astype(str).isin(cache_ids)].reset_index(drop=True)\n    print(f\"Using {len(df)} CTA rows present in cache ({CACHE_DIR})\")\nexcept Exception as e:\n    print(f\"Warning: could not list cache dir {CACHE_DIR}: {e}\")\n\nDEBUG_FAST = False\nif DEBUG_FAST:\n    SUBSET_PER_CLASS = 50\n    df = (df.groupby('Aneurysm Present', group_keys=False)\n            .apply(lambda x: x.sample(min(len(x), SUBSET_PER_CLASS), random_state=42))\n            .reset_index(drop=True))\nelse:\n    df.reset_index(drop=True, inplace=True)\n\nif len(df) == 0:\n    raise ValueError(f\"No CTA rows found in cache after filtering. Check CACHE_DIR ({CACHE_DIR}) and ensure files end with '_img.npy' or '_img.npz'.\")\n\nif USE_CV:\n    df['fold'] = -1\n    skf = StratifiedKFold(n_splits=NUM_FOLDS, shuffle=True, random_state=42)\n    for fold, (_, val_idx) in enumerate(skf.split(df, df['Aneurysm Present'])):\n        df.loc[val_idx, 'fold'] = fold\nelse:\n    df['fold'] = 0\n\n# Training Loop with modifications to handle invalid samples\ncriterion = nn.BCEWithLogitsLoss()\n\ndef train_architecture(arch):\n    for fold in range(NUM_FOLDS if USE_CV else 1):\n        print(f\"🔁 Training {arch} (hybrid) | Fold {fold}\")\n        train_df = df[df['fold'] != fold]\n        val_df   = df[df['fold'] == fold]\n\n        # Version-safe noise: prefer GaussianNoise if available; else fallback to GaussNoise without extra args\n        noise_tf = A.GaussianNoise(var_limit=(1.0, 5.0), mean=0, p=0.2) if hasattr(A, \"GaussianNoise\") else A.GaussNoise(p=0.2)\n        \n        # For hybrid (apply same aug to full and ROI)\n        hybrid_tf = A.Compose(\n            [\n                A.HorizontalFlip(p=0.5),\n                A.Affine(\n                    scale=(0.95, 1.05),\n                    translate_percent={'x': (-0.05, 0.05), 'y': (-0.05, 0.05)},\n                    rotate=(-10, 10),\n                    p=0.5\n                ),\n                noise_tf,\n            ],\n            additional_targets={'image2': 'image'}\n        )\n\n        train_ds = HybridCachedNPYDatasetAlb(train_df, id_col=ID_COL, img_size=IMG_SIZE,\n                                             roi_box_frac=0.15, roi_min_px=24, transform=hybrid_tf)\n        val_ds   = HybridCachedNPYDatasetAlb(val_df,   id_col=ID_COL, img_size=IMG_SIZE,\n                                             roi_box_frac=0.15, roi_min_px=24, transform=None)\n        \n        train_loader = make_loader(train_ds, True,  NUM_WORKERS)\n        val_loader   = make_loader(val_ds,   False, max(1, SAFE_WORKERS//2) if SAFE_WORKERS>0 else 0)\n\n        model = HybridAneurysmModel(base_model_name=arch)\n        if USE_MULTIGPU:\n            model = nn.DataParallel(model)\n        model.to(device)\n\n        optimizer = torch.optim.Adam(model.parameters(), lr=LR)\n        criterion = nn.BCEWithLogitsLoss()\n\n        best_wauc = -1.0\n        for epoch in range(EPOCHS):\n            model.train()\n            running_loss = 0.0\n            pbar = tqdm(train_loader, total=len(train_loader),\n                        desc=f\"Epoch {epoch+1}/{EPOCHS} [train]\", leave=False,\n                       file=sys.stderr)\n            for step, (x_full, x_roi, y, coords) in enumerate(pbar):\n                x_full = x_full.to(device, non_blocking=True)\n                x_roi  = x_roi.to(device, non_blocking=True)\n                y      = y.to(device, non_blocking=True)\n                coords = coords.to(device, non_blocking=True)\n\n                optimizer.zero_grad()\n                logits = model(x_full, x_roi, coords)\n                loss = criterion(logits, y)\n                loss.backward()\n                optimizer.step()\n\n                running_loss += loss.item()\n                if (step + 1) % 10 == 0:\n                    pbar.set_postfix(loss=running_loss / (step + 1))\n\n            model.eval()\n            all_probs, all_targets = [], []\n            for _, row in tqdm(val_df.iterrows(), total=len(val_df), desc=\"val(all-slices)\"):\n                sid = str(row[ID_COL])\n                y   = pd.to_numeric(row[LABEL_COLS], errors=\"coerce\").fillna(0.0).values.astype(np.float32)\n                p   = infer_series_all_slices(model, sid, cache_dir=CACHE_DIR, img_size=IMG_SIZE, offsets=OFFSETS,\n                                            batch_size=16, aggregate=\"max\", use_roi=True)\n                all_probs.append(p)\n                all_targets.append(y)\n\n            all_probs = np.asarray(all_probs, dtype=np.float32)\n            all_targets = np.asarray(all_targets, dtype=np.float32)\n\n            assert all_probs.ndim == 2 and all_targets.ndim == 2\n            assert all_probs.shape[1] == len(LABEL_COLS) and all_targets.shape[1] == len(LABEL_COLS)\n\n            wAUC, ap_auc, others_mean, per_class_auc, skipped = compute_weighted_auc(all_targets, all_probs, LABEL_COLS)\n            print(f\"Skipped {len(skipped)}: {skipped}\")\n            print(f\"Epoch {epoch+1}/{EPOCHS} | train {running_loss/ max(1,len(train_loader)):.4f} | \"\n                  f\"wAUC {wAUC:.4f} | AP {ap_auc:.4f} | others {others_mean:.4f}\")\n            \n            if not np.isnan(wAUC) and wAUC > best_wauc:\n                best_wauc = wAUC\n                torch.save(model.module.state_dict() if USE_MULTIGPU else model.state_dict(),\n                           f\"{arch}_fold{fold}_best_wAUC.pth\")\n\n            del all_probs, all_targets  # and any large tensors like logits/x/y/coords/x_full/x_roi\n            gc.collect()\n            if torch.cuda.is_available():\n                torch.cuda.empty_cache()  # helps fragmentation; doesn't hurt correctness\n\n        model_path = f\"{arch}_hybrid_fold{fold}.pth\"\n        torch.save(model.module.state_dict() if USE_MULTIGPU else model.state_dict(), model_path)\n        print(f\"💾 Saved: {model_path}\")\n\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-12T18:37:44.750315Z","iopub.execute_input":"2025-10-12T18:37:44.75057Z","iopub.status.idle":"2025-10-12T18:37:44.796995Z","shell.execute_reply.started":"2025-10-12T18:37:44.750537Z","shell.execute_reply":"2025-10-12T18:37:44.796387Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Training","metadata":{}},{"cell_type":"code","source":"RUN_TRAIN = True #False to load trained models\n\nif RUN_TRAIN:\n    for arch in ARCHITECTURES:\n        train_architecture(arch)\nelse:\n    models = []\n    for arch in ARCHITECTURES:\n        for fold in range(NUM_FOLDS):\n            p = f\"{arch}_fold{fold}.pth\"\n            if os.path.exists(p):\n                models.append(load_trained_aneurysm_model(p, arch))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-12T18:37:44.797672Z","iopub.execute_input":"2025-10-12T18:37:44.797889Z","execution_failed":"2025-10-13T03:00:07.598Z"}},"outputs":[],"execution_count":null}]}