{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.12.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":117682,"databundleVersionId":15062069,"sourceType":"competition"}],"dockerImageVersionId":31234,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# 🌋 Vesuvius Challenge 📜 Surface Detection\n\n### Official Formula:\n\n    Score = 0.30 × TopoScore + 0.35 × SurfaceDice@τ + 0.35 × VOI_score\n\n### File Format:\n- Input: 3D TIF files (320×320×320 uint8)\n- Labels: 0=background, 1=surface, 2=unlabeled (ignored)\n\n### Usage:\n    \n    gt = load_volume(\"path/to/gt.tif\")\n    pred = load_volume(\"path/to/pred.tif\")\n    result = compute_score(pred, gt)\n    print(f\"Score: {result['score']}\")\n\n### Details\n\nCompetition: https://www.kaggle.com/competitions/vesuvius-challenge-surface-detection\n\nVisualisation at: https://www.kaggle.com/code/jirkaborovec/surface-eda-interactive-img-mask-3d-view\n\nCo-authored by: claude.ai @ Opus 4.5","metadata":{}},{"cell_type":"code","source":"!pip install -q imagecodecs","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:34.194681Z","iopub.execute_input":"2026-01-06T21:10:34.194921Z","iopub.status.idle":"2026-01-06T21:10:40.990771Z","shell.execute_reply.started":"2026-01-06T21:10:34.194897Z","shell.execute_reply":"2026-01-06T21:10:40.989712Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# ⚙️ CONSTANTS\n\n**from competition specification**","metadata":{}},{"cell_type":"code","source":"import time\nimport imagecodecs\nimport tifffile\nimport numpy as np\nfrom scipy import ndimage\nfrom scipy.ndimage import distance_transform_edt, label as cc3d_label\nfrom skimage.measure import euler_number\nfrom typing import Tuple, Dict, Union, List\nfrom pathlib import Path\n\nTAU = 2.0                    # Surface Dice tolerance\nALPHA_VOI = 0.3              # VOI conversion parameter\nBETTI_WEIGHTS = {0: 0.34, 1: 0.33, 2: 0.33}  # TopoScore weights\nW_TOPO = 0.30                # Final score weight for TopoScore\nW_SURFACE_DICE = 0.35        # Final score weight for SurfaceDice\nW_VOI = 0.35                 # Final score weight for VOI\n\nLABEL_BACKGROUND = 0\nLABEL_SURFACE = 1\nLABEL_UNLABELED = 2","metadata":{"trusted":true,"_kg_hide-input":false,"execution":{"iopub.status.busy":"2026-01-06T21:10:40.993066Z","iopub.execute_input":"2026-01-06T21:10:40.993424Z","iopub.status.idle":"2026-01-06T21:10:41.570002Z","shell.execute_reply.started":"2026-01-06T21:10:40.993375Z","shell.execute_reply":"2026-01-06T21:10:41.569045Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# I/O Functions & PreProcessing","metadata":{}},{"cell_type":"code","source":"def load_volume(path: Union[str, Path]) -> np.ndarray:\n    \"\"\"Load 3D volume from TIF file.\"\"\"\n    path = Path(path)\n    if not path.exists():\n        raise FileNotFoundError(f\"Volume not found: {path}\")\n    volume = tifffile.imread(str(path))\n    if volume.ndim != 3:\n        raise ValueError(f\"Expected 3D volume, got shape {volume.shape}\")\n    return volume.astype(np.uint8)\n\n\ndef save_volume(\n    volume: np.ndarray, path: Union[str, Path], compress: bool = True\n):\n    \"\"\"Save 3D volume to TIF file.\"\"\"\n    compression = 'lzw' if compress else None\n    tifffile.imwrite(str(path), volume.astype(np.uint8), compression=compression)\n\n\ndef preprocess(\n    pred: np.ndarray, gt: np.ndarray, \n    ignore_label: int = LABEL_UNLABELED\n) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:\n    \"\"\"Binarize volumes and create valid mask (excluding ignore_label).\"\"\"\n    valid_mask = (gt != ignore_label)\n    pred_bin = ((pred == LABEL_SURFACE) & valid_mask).astype(np.uint8)\n    gt_bin = ((gt == LABEL_SURFACE) & valid_mask).astype(np.uint8)\n    return pred_bin, gt_bin, valid_mask","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.571280Z","iopub.execute_input":"2026-01-06T21:10:41.571884Z","iopub.status.idle":"2026-01-06T21:10:41.580334Z","shell.execute_reply.started":"2026-01-06T21:10:41.571846Z","shell.execute_reply":"2026-01-06T21:10:41.579269Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_components(\n    volume: np.ndarray\n) -> Tuple[np.ndarray, int, Dict[int, int]]:\n    \"\"\"Get connected components with sizes.\"\"\"\n    struct_26 = ndimage.generate_binary_structure(3, 3)\n    labels, n = cc3d_label(volume.astype(bool), structure=struct_26)\n    sizes = {i: int(np.sum(labels == i)) for i in range(1, n + 1)}\n    return labels, n, sizes\n\n\ndef remove_component(\n    volume: np.ndarray, labels: np.ndarray, comp_id: int\n) -> np.ndarray:\n    \"\"\"Remove specific component.\"\"\"\n    result = volume.copy()\n    result[labels == comp_id] = LABEL_BACKGROUND\n    return result\n\n\ndef remove_components(\n    volume: np.ndarray, labels: np.ndarray, comp_ids: List[int]\n) -> np.ndarray:\n    \"\"\"Remove multiple components.\"\"\"\n    result = volume.copy()\n    for cid in comp_ids:\n        result[labels == cid] = LABEL_BACKGROUND\n    return result\n\n\ndef add_hole_in_component(\n    volume: np.ndarray,\n    labels: np.ndarray,\n    comp_id: int,\n    hole_size: int = 25\n) -> Tuple[np.ndarray, int]:\n    \"\"\"\n    Add a rectangular hole in the center of a specific component.\n    \n    Creates a topology break (adds β₁ tunnel) without affecting other components.\n    \n    Args:\n        volume: Label volume\n        labels: Connected component labels\n        comp_id: Which component to add hole to\n        hole_size: Half-size of hole in Y and Z dimensions\n        \n    Returns:\n        (modified_volume, voxels_removed)\n    \"\"\"\n    comp_mask = (labels == comp_id)\n    coords = np.where(comp_mask)\n    \n    if len(coords[0]) == 0:\n        return volume.copy(), 0\n    \n    # Find center of component\n    z_center = int(np.mean(coords[0]))\n    y_center = int(np.mean(coords[1]))\n    \n    # Create rectangular hole region\n    D, H, W = volume.shape\n    z_grid = np.arange(D)[:, None, None]\n    y_grid = np.arange(H)[None, :, None]\n    \n    hole_region = (\n        (z_grid >= z_center - hole_size) & \n        (z_grid <= z_center + hole_size) &\n        (y_grid >= y_center - hole_size) & \n        (y_grid <= y_center + hole_size)\n    )\n    \n    # Only remove voxels in BOTH hole region AND this component\n    hole_mask = hole_region & comp_mask\n    \n    result = volume.copy()\n    result[hole_mask] = LABEL_BACKGROUND\n    \n    return result, int(np.sum(hole_mask))\n\n\ndef analyze_volume(volume: np.ndarray, name: str = \"Volume\") -> Dict:\n    \"\"\"Analyze label volume and return statistics.\"\"\"\n    surface = (volume == LABEL_SURFACE).astype(np.uint8)\n    labels, n_comp, sizes = get_components(surface)\n    betti = compute_betti_numbers(surface)\n    \n    return {\n        'shape': volume.shape,\n        'n_surface_voxels': int(np.sum(surface)),\n        'n_components': n_comp,\n        'betti': betti,\n        'component_labels': labels,\n        'component_sizes': sizes\n    }\n\n\ndef add_hole_in_component(\n    volume: np.ndarray,\n    labels: np.ndarray,\n    comp_id: int,\n    hole_size: int = 25,\n    label_bg: int = 0,\n) -> Tuple[np.ndarray, int]:\n    \"\"\"Add a rectangular hole in the center of a specific component.\n    \n    Creates a topology break (adds β₁ tunnel) without affecting other components.\n    \n    Args:\n        volume: Label volume\n        labels: Connected component labels\n        comp_id: Which component to add hole to\n        hole_size: Half-size of hole in Y and Z dimensions\n        \n    Returns:\n        (modified_volume, voxels_removed)\n    \"\"\"\n    comp_mask = (labels == comp_id)\n    coords = np.where(comp_mask)\n    \n    if len(coords[0]) == 0:\n        return volume.copy(), 0\n    \n    # Find center of component\n    z_center = int(np.mean(coords[0]))\n    y_center = int(np.mean(coords[1]))\n    \n    # Create rectangular hole region\n    D, H, W = volume.shape\n    z_grid = np.arange(D)[:, None, None]\n    y_grid = np.arange(H)[None, :, None]\n    \n    hole_region = (\n        (z_grid >= z_center - hole_size) & \n        (z_grid <= z_center + hole_size) &\n        (y_grid >= y_center - hole_size) & \n        (y_grid <= y_center + hole_size)\n    )\n    \n    # Only remove voxels in BOTH hole region AND this component\n    hole_mask = hole_region & comp_mask\n    \n    result = volume.copy()\n    result[hole_mask] = label_bg\n    \n    return result, int(np.sum(hole_mask))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.581705Z","iopub.execute_input":"2026-01-06T21:10:41.582064Z","iopub.status.idle":"2026-01-06T21:10:41.601715Z","shell.execute_reply.started":"2026-01-06T21:10:41.582029Z","shell.execute_reply":"2026-01-06T21:10:41.600733Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 1 🎲 Surface Dice @ τ","metadata":{}},{"cell_type":"code","source":"def extract_surface(volume: np.ndarray) -> np.ndarray:\n    \"\"\"Extract boundary voxels using 6-connectivity erosion.\"\"\"\n    if np.sum(volume) == 0:\n        return np.zeros_like(volume, dtype=np.uint8)\n    struct_6 = ndimage.generate_binary_structure(3, 1)\n    eroded = ndimage.binary_erosion(volume.astype(bool), structure=struct_6)\n    return volume.astype(np.uint8) - eroded.astype(np.uint8)\n\n\ndef compute_surface_dice(\n    pred: np.ndarray, gt: np.ndarray, tau: float = TAU\n) -> float:\n    \"\"\"\n    Surface Dice @ τ: fraction of surface points within tolerance.\n    \n    Edge cases: both empty → 1.0, one empty → 0.0\n    \"\"\"\n    pred_empty = np.sum(pred) == 0\n    gt_empty = np.sum(gt) == 0\n    \n    if pred_empty and gt_empty:\n        return 1.0\n    if pred_empty or gt_empty:\n        return 0.0\n    \n    pred_surface = extract_surface(pred)\n    gt_surface = extract_surface(gt)\n    \n    n_pred = np.sum(pred_surface)\n    n_gt = np.sum(gt_surface)\n    \n    if n_pred == 0 or n_gt == 0:\n        return 0.0\n    \n    # Distance from pred surface to GT SURFACE (not volume!)\n    dist_to_gt_surface = distance_transform_edt(~gt_surface.astype(bool))\n    pred_matched = np.sum(dist_to_gt_surface[pred_surface > 0] <= tau)\n    \n    # Distance from GT surface to pred SURFACE (not volume!)\n    dist_to_pred_surface = distance_transform_edt(~pred_surface.astype(bool))\n    gt_matched = np.sum(dist_to_pred_surface[gt_surface > 0] <= tau)\n    \n    return float((pred_matched + gt_matched) / (n_pred + n_gt))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.604528Z","iopub.execute_input":"2026-01-06T21:10:41.605095Z","iopub.status.idle":"2026-01-06T21:10:41.621242Z","shell.execute_reply.started":"2026-01-06T21:10:41.605049Z","shell.execute_reply":"2026-01-06T21:10:41.620217Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 2 📦 VOI Score","metadata":{}},{"cell_type":"code","source":"def compute_voi(\n    pred_labels: np.ndarray, gt_labels: np.ndarray\n) -> Tuple[float, float, float]:\n    \"\"\"Compute Variation of Information: (VOI_split, VOI_merge, VOI_total).\"\"\"\n    n = pred_labels.size\n    if n == 0:\n        return 0.0, 0.0, 0.0\n    \n    max_pred = int(np.max(pred_labels)) + 1\n    max_gt = int(np.max(gt_labels)) + 1\n    \n    contingency = np.zeros((max_pred, max_gt), dtype=np.float64)\n    for p, g in zip(pred_labels.ravel(), gt_labels.ravel()):\n        contingency[p, g] += 1\n    \n    p_ij = contingency / n\n    p_i = np.sum(p_ij, axis=1)\n    p_j = np.sum(p_ij, axis=0)\n    \n    voi_split = 0.0\n    voi_merge = 0.0\n    \n    for i in range(max_pred):\n        for j in range(max_gt):\n            if p_ij[i, j] > 0:\n                if p_i[i] > 0:\n                    voi_split -= p_ij[i, j] * np.log2(p_ij[i, j] / p_i[i])\n                if p_j[j] > 0:\n                    voi_merge -= p_ij[i, j] * np.log2(p_ij[i, j] / p_j[j])\n    \n    return float(voi_split), float(voi_merge), float(voi_split + voi_merge)\n\n\ndef compute_voi_score(pred: np.ndarray, gt: np.ndarray) -> float:\n    \"\"\"VOI Score = 1 / (1 + α × VOI_total). Uses 26-connectivity.\"\"\"\n    if np.sum(pred) == 0 and np.sum(gt) == 0:\n        return 1.0\n    if np.sum(pred) == 0 or np.sum(gt) == 0:\n        return 0.0\n    \n    struct_26 = ndimage.generate_binary_structure(3, 3)\n    pred_labels, _ = cc3d_label(pred.astype(bool), structure=struct_26)\n    gt_labels, _ = cc3d_label(gt.astype(bool), structure=struct_26)\n    \n    union_fg = (pred.astype(bool) | gt.astype(bool))\n    _, _, voi_total = compute_voi(pred_labels[union_fg], gt_labels[union_fg])\n    \n    return float(1.0 / (1.0 + ALPHA_VOI * voi_total))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.622664Z","iopub.execute_input":"2026-01-06T21:10:41.623069Z","iopub.status.idle":"2026-01-06T21:10:41.638062Z","shell.execute_reply.started":"2026-01-06T21:10:41.623027Z","shell.execute_reply":"2026-01-06T21:10:41.637107Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 3 🏔️ TOPO Score","metadata":{}},{"cell_type":"code","source":"def compute_betti_numbers(volume: np.ndarray) -> Tuple[int, int, int]:\n    \"\"\"Compute Betti numbers (β₀, β₁, β₂) for 3D binary volume.\"\"\"\n    if np.sum(volume) == 0:\n        return (0, 0, 0)\n    \n    vol_bool = volume.astype(bool)\n    struct_26 = ndimage.generate_binary_structure(3, 3)\n    _, beta0 = cc3d_label(vol_bool, structure=struct_26)\n    \n    # Cavities (enclosed background components)\n    struct_6 = ndimage.generate_binary_structure(3, 1)\n    bg_labels, n_bg = cc3d_label(~vol_bool, structure=struct_6)\n    \n    boundary_labels = set()\n    for face in [bg_labels[0,:,:], bg_labels[-1,:,:],\n                 bg_labels[:,0,:], bg_labels[:,-1,:],\n                 bg_labels[:,:,0], bg_labels[:,:,-1]]:\n        boundary_labels.update(face.ravel())\n    boundary_labels.discard(0)\n    beta2 = n_bg - len(boundary_labels)\n    \n    # β₁ from Euler characteristic\n    chi = euler_number(vol_bool, connectivity=3)\n    beta1 = max(0, beta0 - chi + beta2)\n    \n    return (int(beta0), int(beta1), int(beta2))\n\n\ndef compute_betti_f1(pred_betti: int, gt_betti: int) -> float:\n    \"\"\"F1 score for Betti number matching.\"\"\"\n    if pred_betti == 0 and gt_betti == 0:\n        return 1.0\n    if pred_betti == 0 or gt_betti == 0:\n        return 0.0\n    matched = min(pred_betti, gt_betti)\n    precision = matched / pred_betti\n    recall = matched / gt_betti\n    if precision + recall == 0:\n        return 0.0\n    return 2.0 * precision * recall / (precision + recall)\n\n\ndef compute_topo_score(\n    pred: np.ndarray, gt: np.ndarray\n) -> Tuple[float, Tuple, Tuple]:\n    \"\"\"TopoScore: weighted Betti F1. Returns (score, pred_betti, gt_betti).\"\"\"\n    pred_betti = compute_betti_numbers(pred)\n    gt_betti = compute_betti_numbers(gt)\n    \n    f1_scores = {}\n    active_weights = {}\n    \n    for k in range(3):\n        pb, gb = pred_betti[k], gt_betti[k]\n        if pb > 0 or gb > 0:\n            f1_scores[k] = compute_betti_f1(pb, gb)\n            active_weights[k] = BETTI_WEIGHTS[k]\n    \n    if not active_weights:\n        return 1.0, pred_betti, gt_betti\n    \n    total_weight = sum(active_weights.values())\n    score = sum(f1_scores[k] * active_weights[k] for k in active_weights) / total_weight\n    \n    return float(score), pred_betti, gt_betti","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.639873Z","iopub.execute_input":"2026-01-06T21:10:41.640240Z","iopub.status.idle":"2026-01-06T21:10:41.659122Z","shell.execute_reply.started":"2026-01-06T21:10:41.640203Z","shell.execute_reply":"2026-01-06T21:10:41.658023Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# ⚖️ MAIN 🏆 Scoring function","metadata":{}},{"cell_type":"code","source":"from functools import wraps\n\ndef time_log(func):\n    @wraps(func)\n    def wrapper(*args, **kwargs):\n        start = time.perf_counter()\n        result = func(*args, **kwargs)\n        end = time.perf_counter()\n        print(f\"[Timer] {func.__name__:20} | {end - start:.4f}s\")\n        return result\n    return wrapper","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.660855Z","iopub.execute_input":"2026-01-06T21:10:41.661241Z","iopub.status.idle":"2026-01-06T21:10:41.679023Z","shell.execute_reply.started":"2026-01-06T21:10:41.661205Z","shell.execute_reply":"2026-01-06T21:10:41.677992Z"},"_kg_hide-input":true,"_kg_hide-output":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def compute_score(\n    pred: np.ndarray,\n    gt: np.ndarray,\n    tau: float = TAU,\n    ignore_label: int = LABEL_UNLABELED,\n    verbose: bool = False\n) -> Dict:\n    \"\"\"Compute official Vesuvius Surface Detection competition score.\n    \n    Score = 0.30 × TopoScore + 0.35 × SurfaceDice@τ + 0.35 × VOI_score\n    \"\"\"\n    # 1. Preprocess\n    t0 = time.perf_counter()\n    pred_bin, gt_bin, _ = preprocess(pred, gt, ignore_label)\n    t_pre = time.perf_counter() - t0\n\n    # 2. Surface Dice\n    t1 = time.perf_counter()\n    surface_dice = compute_surface_dice(pred_bin, gt_bin, tau)\n    t_dice = time.perf_counter() - t1\n\n    # 3. VOI\n    t2 = time.perf_counter()\n    voi_score = compute_voi_score(pred_bin, gt_bin)\n    t_voi = time.perf_counter() - t2\n\n    # 4. Topo Score\n    t3 = time.perf_counter()\n    topo_score, pred_betti, gt_betti = compute_topo_score(pred_bin, gt_bin)\n    t_topo = time.perf_counter() - t3\n    \n    score = W_TOPO * topo_score + W_SURFACE_DICE * surface_dice + W_VOI * voi_score\n    \n    result = {\n        'score': round(score, 6),\n        'surface_dice': round(surface_dice, 6),\n        'voi_score': round(voi_score, 6),\n        'topo_score': round(topo_score, 6),\n        'pred_betti': pred_betti,\n        'gt_betti': gt_betti\n    }\n    \n    if verbose:\n        print(f\"\\n--- Timing Breakdown ---\")\n        print(f\"Preprocess:   {t_pre:.4f}s\")\n        print(f\"Surface Dice: {t_dice:.4f}s\")\n        print(f\"VOI Score:    {t_voi:.4f}s\")\n        print(f\"Topo Score:   {t_topo:.4f}s\")\n        print(f\"--- Particular Metrics ---\")\n        print(f\"SurfaceDice@{tau}: {surface_dice:.6f}\")\n        print(f\"VOI_score:        {voi_score:.6f}\")\n        print(f\"TopoScore:        {topo_score:.6f}\")\n        print(f\"  Pred β: {pred_betti}, GT β: {gt_betti}\")\n        print(f\">>> SCORE:        {score:.6f}\\n\")\n    \n    return result\n\n\ndef score_from_files(\n    pred_path: Union[str, Path],\n    gt_path: Union[str, Path],\n    tau: float = TAU,\n    verbose: bool = True\n) -> Dict:\n    \"\"\"Compute score directly from TIF file paths.\"\"\"\n    if verbose:\n        print(f\"Loading: {pred_path}\")\n    pred = load_volume(pred_path)\n    \n    if verbose:\n        print(f\"Loading: {gt_path}\")\n    gt = load_volume(gt_path)\n    \n    if pred.shape != gt.shape:\n        raise ValueError(f\"Shape mismatch: pred {pred.shape} vs gt {gt.shape}\")\n    \n    if verbose:\n        print(f\"Shape: {gt.shape}\")\n        print(\"-\" * 50)\n    \n    return compute_score(pred, gt, tau=tau, verbose=verbose)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.680191Z","iopub.execute_input":"2026-01-06T21:10:41.680567Z","iopub.status.idle":"2026-01-06T21:10:41.702752Z","shell.execute_reply.started":"2026-01-06T21:10:41.680540Z","shell.execute_reply":"2026-01-06T21:10:41.701773Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# 🚀 Real 📽️ DEMO","metadata":{}},{"cell_type":"code","source":"print(\"=\" * 70)\nprint(\"VESUVIUS SURFACE DETECTION - KAGGLE COMPETITION METRIC\")\nprint(\"=\" * 70)\nprint()\nprint(\"Score = 0.30×TopoScore + 0.35×SurfaceDice@2.0 + 0.35×VOI_score\")\nprint()\n\n# # Score from files\n# result = score_from_files(sys.argv[1], sys.argv[2], verbose=True)\n\n# Analyze single file and run instance tests\ngt_path = \"/kaggle/input/vesuvius-challenge-surface-detection/train_labels/602831951.tif\"\nprint(f\"Loading: {gt_path}\")\ngt = load_volume(gt_path)\n\n# Analyze\ninfo = analyze_volume(gt)\nprint(f\"\\nVolume: {info['shape']}\")\nprint(f\"Surface voxels: {info['n_surface_voxels']:,}\")\nprint(f\"Components: {info['n_components']}\")\nprint(f\"Betti: β₀={info['betti'][0]}, β₁={info['betti'][1]}, β₂={info['betti'][2]}\")\n\n# Sort components by size\nsorted_comps = sorted(info['component_sizes'].items(), key=lambda x: -x[1])\nprint(f\"\\nTop 5 components:\")\nfor i, (cid, size) in enumerate(sorted_comps[:5]):\n    pct = 100 * size / info['n_surface_voxels']\n    print(f\"  #{i+1} (ID {cid}): {size:,} voxels ({pct:.1f}%)\")\n\n# Run test cases\nprint(\"\\n\" + \"=\" * 70)\nprint(\"INSTANCE-BASED TEST CASES\")\nprint(\"=\" * 70)\n\nlabels = info['component_labels']\nlargest_id = sorted_comps[0][0]\n\ntests = [\n    (\"Perfect Match\", gt.copy()),\n    (\"Remove 1 Surface\", remove_component(gt.copy(), labels, sorted_comps[0][0])),\n    (\"Remove 3 Surfaces\", remove_components(gt.copy(), labels, [c[0] for c in sorted_comps[:3]])),\n    # (\"Remove Half\", remove_components(gt.copy(), labels, [c[0] for c in sorted_comps[:len(sorted_comps)//2]])),\n    (\"Empty Prediction\", np.zeros_like(gt)),\n]\nfor hole in (25, 50, 100, 150):\n    # Create hole in largest component\n    pred_with_hole, voxels_removed = add_hole_in_component(gt.copy(), labels, largest_id, hole_size=hole)\n    tests.append((f\"Hole in Surface (−{voxels_removed:,} vox)\", pred_with_hole)) \n\nresults = []\nfor name, pred in tests:\n    r = compute_score(pred, gt, verbose=True)\n    results.append(f\"{name:<35}\" \\\n        f\" {r['surface_dice']:>10.4f} {r['voi_score']:>10.4f} {r['topo_score']:>10.4f} {r['score']:>10.4f}\")\n\nprint(f\"\\n{'Case':<35} {'SurfDice':>10} {'VOI':>10} {'Topo':>10} {'SCORE':>10}\")\nprint(\"-\" * 70)\nfor res in results:\n    print(res)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2026-01-06T21:10:41.703935Z","iopub.execute_input":"2026-01-06T21:10:41.704226Z","iopub.status.idle":"2026-01-06T21:13:32.590333Z","shell.execute_reply.started":"2026-01-06T21:10:41.704203Z","shell.execute_reply":"2026-01-06T21:13:32.589299Z"}},"outputs":[],"execution_count":null}]}