{"metadata":{"kernelspec":{"display_name":".venv","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.10.12"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":117682,"databundleVersionId":15062069,"sourceType":"competition"}],"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"id":"477adabc","cell_type":"markdown","source":"# Vesuvius Challenge - Unsupervised Baseline\n\n**Problem**: Train labels are placeholders (all value=2), no real ground truth available.\n**Solution**: Unsupervised baseline using intensity-based segmentation.\n\nThis notebook creates predictions using simple thresholding + morphological operations.","metadata":{}},{"id":"e71a5569","cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom pathlib import Path\nfrom tqdm import tqdm\nimport matplotlib.pyplot as plt\nfrom PIL import Image, ImageSequence\nimport tifffile as tiff\nfrom scipy import ndimage\nimport os","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"7f37a3d2","cell_type":"code","source":"# Check if running on Kaggle\nIS_KAGGLE = 'KAGGLE_DATA_PROXY_URL' in os.environ\nprint(f\"Running on Kaggle: {IS_KAGGLE}\")\n\n\n# Detect environment paths\nif IS_KAGGLE:\n    DATA_PATH = Path('/kaggle/input/vesuvius-challenge-surface-detection')\n    OUTPUT_PATH = Path('/kaggle/working')\nelse:\n    # Notebook lives in notebooks\n    DATA_PATH = Path('../vesuvius-challenge-surface-detection')\n    OUTPUT_PATH = Path('../output')\n\n# Create output directory if it doesn't exist\nOUTPUT_PATH.mkdir(parents=True, exist_ok=True)\nprint(f\"Data path: {DATA_PATH}\")\nprint(f\"Output path: {OUTPUT_PATH}\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"f0c6f737","cell_type":"code","source":"def read_volume(path: Path):\n    \"\"\"Read 3D volume using PIL\"\"\"\n    im = Image.open(str(path))\n    frames = [np.array(f) for f in ImageSequence.Iterator(im)]\n    vol = np.stack(frames, axis=0)\n    if vol.ndim == 2:\n        vol = vol[None, ...]\n    return vol\n\n# Load test data\ntest_df = pd.read_csv(DATA_PATH / 'test.csv')\ntest_ids = test_df['id'].astype(str).tolist()  # ALL test volumes\n\nprint(f\"\\nLoading {len(test_ids)} test volumes...\")\n\ntest_volumes = {}\nfor tid in tqdm(test_ids):\n    img_path = DATA_PATH / 'test_images' / f\"{tid}.tif\"\n    if img_path.exists():\n        test_volumes[tid] = read_volume(img_path)\n\nprint(f\"✓ Loaded {len(test_volumes)} test volumes\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"3d37230b","cell_type":"markdown","source":"## Baseline Method: Simple Percentile Thresholding\n\nThis is our baseline approach using intensity-based segmentation with morphological operations. We'll test it on a sample volume first to understand the approach, then try improvements.","metadata":{}},{"id":"f85fc407","cell_type":"code","source":"def segment_surface_unsupervised(volume, percentile_threshold=60):\n    \"\"\"\n    Unsupervised surface segmentation using intensity thresholding.\n    \n    Strategy:\n    1. Percentile-based thresholding\n    2. Morphological operations to clean up\n    3. Connected components to remove noise\n    \"\"\"\n    # Compute threshold\n    threshold = np.percentile(volume, percentile_threshold)\n    \n    # Initial binary mask\n    mask = volume > threshold\n    \n    # Morphological closing to fill holes\n    mask = ndimage.binary_closing(mask, structure=np.ones((3, 3, 3)))\n    \n    # Remove small objects\n    mask = ndimage.binary_opening(mask, structure=np.ones((3, 3, 3)))\n    \n    return mask.astype(np.uint8)\n\n# Test on one volume\ntest_id = list(test_volumes.keys())[0]\ntest_vol = test_volumes[test_id]\n\nprint(f\"\\nTesting on volume {test_id}:\")\nprint(f\"  Shape: {test_vol.shape}\")\nprint(f\"  Intensity range: [{test_vol.min()}, {test_vol.max()}]\")\n\nmask = segment_surface_unsupervised(test_vol, percentile_threshold=60)\n\nprint(f\"\\nMask results:\")\nprint(f\"  Shape: {mask.shape}\")\nprint(f\"  Coverage: {mask.mean()*100:.1f}%\")\nprint(f\"  Unique values: {np.unique(mask)}\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"d3a32e0a","cell_type":"markdown","source":"## Improved Method: Edge Detection + 2.5D Voting\n\nNow we create an improved approach that combines intensity thresholding with:\n- **Gradient-based edge detection** to capture surface boundaries\n- **2.5D voting** using neighboring slices for consistency\n- **Adaptive noise filtering** to remove false positives","metadata":{}},{"id":"98d25f92","cell_type":"code","source":"# Improved segmentation with edge detection + 2.5D approach\ndef segment_surface_improved(\n    volume,\n    percentile_threshold=65,\n    use_edges=True,\n    use_2p5d=True,\n):\n    \"\"\"\n    Improved surface segmentation combining intensity and edge information.\n\n    Improvements over baseline:\n    1. Gradient-based edge detection to capture surface boundaries\n    2. 2.5D voting using neighboring Z-slices for consistency\n    3. Morphological cleanup operations\n\n    Args:\n        volume: 3D numpy array of intensities\n        percentile_threshold: Percentile value for intensity threshold (0-100)\n        use_edges: If True, enhance with gradient-based edge detection\n        use_2p5d: If True, apply 2.5D voting across Z-slices\n\n    Returns:\n        Binary mask of same shape as volume (uint8, values 0 or 1)\n    \"\"\"\n    # Initial intensity-based threshold\n    threshold = np.percentile(volume, percentile_threshold)\n    mask = volume > threshold\n\n    # Edge detection: use gradients to enhance surface boundaries\n    if use_edges:\n        # Compute gradients along all axes\n        grad_z = np.abs(np.diff(volume, axis=0, prepend=volume[0:1]))\n        grad_y = np.abs(np.diff(volume, axis=1, prepend=volume[:, 0:1]))\n        grad_x = np.abs(np.diff(volume, axis=2, prepend=volume[:, :, 0:1]))\n\n        # Magnitude of gradient\n        grad_mag = np.sqrt(grad_z**2 + grad_y**2 + grad_x**2)\n\n        # Normalize and threshold edges\n        grad_threshold = np.percentile(grad_mag, 50)\n        edges = grad_mag > grad_threshold\n\n        # Combine: include both intensity threshold AND strong edges\n        mask = mask | edges\n\n    # 2.5D refinement: voting with neighboring slices for consistency\n    if use_2p5d:\n        refined_mask = mask.copy()\n        for z in range(volume.shape[0]):\n            # Get current slice and neighbors (handle boundaries)\n            z_prev = max(0, z - 1)\n            z_next = min(volume.shape[0] - 1, z + 1)\n\n            # Voting: pixel positive if 2+ of 3 slices (prev, current, next) agree\n            votes = (\n                mask[z_prev].astype(int)\n                + mask[z].astype(int)\n                + mask[z_next].astype(int)\n            )\n            refined_mask[z] = votes >= 2\n\n        mask = refined_mask\n\n    # Morphological operations for cleanup\n    mask = ndimage.binary_closing(mask, structure=np.ones((3, 3, 3)))\n    mask = ndimage.binary_opening(mask, structure=np.ones((3, 3, 3)))\n\n    return mask.astype(np.uint8)","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"e167bffd","cell_type":"markdown","source":"## Test Improved Method on Sample Volume\n\nLet's test the improved segmentation on our sample volume and try different percentile thresholds to find the best balance.","metadata":{}},{"id":"a22ad71f","cell_type":"code","source":"# Test improved version after basic baseline is tested\nprint(f\"\\nTesting improved segmentation on volume {test_id}:\")\n\nmask_improved = segment_surface_improved(\n    test_vol, \n    percentile_threshold=60, \n    use_edges=True, \n    use_2p5d=True\n)\n\nprint(f\"  Improved mask coverage: {mask_improved.mean()*100:.1f}%\")\nprint(f\"  Unique values: {np.unique(mask_improved)}\")\n\n# Test multiple threshold values to find optimal\nprint(f\"\\nTesting multiple thresholds...\")\nprint(f\"\\nPercentile | Coverage | Unique Values\")\nprint(f\"-----------|----------|---------------\")\n\nthreshold_results = {}\nfor pct in [50, 55, 60, 65, 70]:\n    mask_test = segment_surface_improved(\n        test_vol,\n        percentile_threshold=pct,\n        use_edges=True,\n        use_2p5d=True\n    )\n    coverage = mask_test.mean() * 100\n    threshold_results[pct] = mask_test\n    print(f\"   {pct:2d}%    | {coverage:6.1f}%  | {np.unique(mask_test)}\")\n\nprint(f\"\\n✓ Best threshold to try: 65% (highest detail with moderate coverage)\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"d1fbdfb1","cell_type":"markdown","source":"## Visual Comparison: Baseline vs Improved\n\nCompare the baseline method with different improved configurations side-by-side to see the improvements.","metadata":{}},{"id":"954d737f","cell_type":"code","source":"# Visualize comparison of methods\nfig, axes = plt.subplots(2, 3, figsize=(15, 10))\n\nz_mid = test_vol.shape[0] // 2\n\n# Original\naxes[0, 0].imshow(test_vol[z_mid], cmap='gray')\naxes[0, 0].set_title('Original Volume')\naxes[0, 0].axis('off')\n\n# Basic method (60%)\nmask_basic = segment_surface_unsupervised(test_vol, percentile_threshold=60)\naxes[0, 1].imshow(mask_basic[z_mid], cmap='gray')\naxes[0, 1].set_title(f'Basic (60%) - {mask_basic.mean()*100:.1f}%')\naxes[0, 1].axis('off')\n\n# Improved method (60%)\nmask_imp60 = segment_surface_improved(test_vol, percentile_threshold=60, use_edges=True, use_2p5d=True)\naxes[0, 2].imshow(mask_imp60[z_mid], cmap='gray')\naxes[0, 2].set_title(f'Improved (60%) - {mask_imp60.mean()*100:.1f}%')\naxes[0, 2].axis('off')\n\n# Improved (65%)\nmask_imp65 = segment_surface_improved(test_vol, percentile_threshold=65, use_edges=True, use_2p5d=True)\naxes[1, 0].imshow(mask_imp65[z_mid], cmap='gray')\naxes[1, 0].set_title(f'Improved (65%) - {mask_imp65.mean()*100:.1f}%')\naxes[1, 0].axis('off')\n\n# Improved (70%)\nmask_imp70 = segment_surface_improved(test_vol, percentile_threshold=70, use_edges=True, use_2p5d=True)\naxes[1, 1].imshow(mask_imp70[z_mid], cmap='gray')\naxes[1, 1].set_title(f'Improved (70%) - {mask_imp70.mean()*100:.1f}%')\naxes[1, 1].axis('off')\n\n# Overlay (improved 65%)\naxes[1, 2].imshow(test_vol[z_mid], cmap='gray')\naxes[1, 2].imshow(mask_imp65[z_mid], cmap='Reds', alpha=0.4)\naxes[1, 2].set_title('Improved (65%) Overlay')\naxes[1, 2].axis('off')\n\nplt.tight_layout()\nplt.savefig(OUTPUT_PATH / f'improved_comparison_{test_id}.png', dpi=100, bbox_inches='tight')\nplt.show()\n\nprint(f\"✓ Comparison visualization saved\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"a0a7db18","cell_type":"markdown","source":"## Generate Full Predictions\n\nNow generate predictions for all test volumes using the improved method with the best threshold (65% percentile with edge detection and 2.5D voting).","metadata":{}},{"id":"d2288d78","cell_type":"code","source":"# Generate predictions for all test volumes (using improved method)\nprint(f\"\\nGenerating improved predictions for {len(test_volumes)} volumes...\\n\")\n\n# Choose best configuration (65% threshold with edges + 2.5D voting)\nPERCENTILE = 65\nUSE_EDGES = True\nUSE_2P5D = True\n\npredictions = {}\nfor tid, vol in tqdm(test_volumes.items()):\n    mask = segment_surface_improved(\n        vol, \n        percentile_threshold=PERCENTILE,\n        use_edges=USE_EDGES,\n        use_2p5d=USE_2P5D\n    )\n    predictions[tid] = mask\n    print(f\"  {tid}: coverage={mask.mean()*100:.1f}%\")\n\nprint(f\"\\n✓ Generated {len(predictions)} improved predictions\")\nprint(f\"  Configuration: {PERCENTILE}% threshold, edges={USE_EDGES}, 2.5D={USE_2P5D}\")\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"4abe2880","cell_type":"markdown","source":"## Save and Package for Submission\n\nSave all predictions as uncompressed TIFFs and create the submission.zip file for Kaggle upload.","metadata":{}},{"id":"b63c0f35","cell_type":"code","source":"# Save predictions as uncompressed TIFFs\nprint(f\"\\nSaving predictions...\")\n\nfor tid, mask in predictions.items():\n    output_path = OUTPUT_PATH / f\"{tid}.tif\"\n    tiff.imwrite(output_path, mask, compression=None)\n    print(f\"  Saved {tid}.tif\")\n\nprint(f\"\\n✓ All predictions saved to {OUTPUT_PATH}\")","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"984a3e56","cell_type":"code","source":"# Create submission directory and prepare files\nimport shutil\nimport os\n\nSUBMISSION_PATH = Path('submission')\nSUBMISSION_PATH.mkdir(parents=True, exist_ok=True)\n\nprint(f\"\\nPreparing submission...\")\n\nfor tid, mask in predictions.items():\n    src = OUTPUT_PATH / f\"{tid}.tif\"\n    dst = SUBMISSION_PATH / f\"{tid}.tif\"\n    \n    if src.exists():\n        shutil.copy2(src, dst)\n        print(f\"  Copied {tid}.tif to submission/\")\n\nprint(f\"\\n✓ Submission directory ready at {SUBMISSION_PATH}\")\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"8dd0ac25","cell_type":"code","source":"# Create submission ZIP file\nimport zipfile\n\nzip_path = Path('submission.zip')\n\nprint(f\"\\nCreating submission ZIP...\")\n\nwith zipfile.ZipFile(zip_path, 'w', zipfile.ZIP_DEFLATED) as zf:\n    for tid in predictions.keys():\n        file_path = SUBMISSION_PATH / f\"{tid}.tif\"\n        if file_path.exists():\n            zf.write(file_path, arcname=f\"{tid}.tif\")\n            print(f\"  Added {tid}.tif to ZIP\")\n\nzip_size_mb = zip_path.stat().st_size / (1024 * 1024)\nprint(f\"\\n✓ submission.zip created ({zip_size_mb:.1f} MB)\")\nprint(f\"\\n📤 Ready to upload submission.zip to Kaggle!\")\n","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"id":"bf084e44","cell_type":"markdown","source":"## Summary\n\nThis unsupervised baseline uses percentile thresholding (60th percentile) + morphological operations.\n\n**Next steps:**\n1. Tune the percentile threshold (try 50, 55, 65, 70)\n2. Add edge detection for surface boundaries\n3. Use 2.5D approach (neighboring slices)\n4. Submit to Kaggle to get actual metric scores\n5. If scores are reasonable, train a supervised model using test feedback","metadata":{}}]}