{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":99552,"databundleVersionId":13851420,"isSourceIdPinned":false,"sourceType":"competition"},{"sourceId":592755,"sourceType":"modelInstanceVersion","isSourceIdPinned":false,"modelInstanceId":443543,"modelId":460052}],"dockerImageVersionId":31090,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"# Install ultralytics (YOLOv8), pydicom, and OpenCV\n!pip install --upgrade pip\n!pip install ultralytics pydicom opencv-python-headless","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-09-28T23:42:17.36324Z","iopub.execute_input":"2025-09-28T23:42:17.363884Z","iopub.status.idle":"2025-09-28T23:49:01.949661Z","shell.execute_reply.started":"2025-09-28T23:42:17.363853Z","shell.execute_reply":"2025-09-28T23:49:01.948892Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"from glob import glob\nimport pydicom as dicom #for dicom files\nimport nibabel as nib #for nii files\n\n\nimport os\nimport shutil\nimport gc\nfrom collections import defaultdict\nfrom typing import Tuple, List\n\n\nimport numpy as np\nimport pandas as pd\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nimport polars as pl\nimport pydicom\nfrom scipy import ndimage\nfrom sklearn.preprocessing import StandardScaler\n\n\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.data import Dataset, DataLoader\nimport torch.optim as optim\n\n\nimport kaggle_evaluation.rsna_inference_server\n\n\nimport warnings\nwarnings.filterwarnings('ignore')\n\ntrain_images = glob(\"/kaggle/input/rsna-intracranial-aneurysm-detection/series/1.2.826.0.1.3680043.8.498.10004044428023505108375152878107656647/*\")\n\npath = '/kaggle/input/rsna-intracranial-aneurysm-detection/segmentations/1.2.826.0.1.3680043.8.498.10035643165968342618460849823699311381.nii'\n\ntrain_df = pd.read_csv(\"/kaggle/input/rsna-intracranial-aneurysm-detection/train.csv\")\nlabel_df = pd.read_csv(\"/kaggle/input/rsna-intracranial-aneurysm-detection/train_localizers.csv\")\n\nplt.style.use('default')\nfig, axes = plt.subplots(4,4, figsize=(12,12))\ntrain_images\nfor i, ax in enumerate(axes.reshape(-1)):\n    img_path = train_images[i]\n    img = dicom.dcmread(img_path)  \n    ax.imshow(img.pixel_array)\nplt.show()\n\nimg = nib.load(path).get_fdata()\nimg.shape\n\n\nplt.style.use('default')\nfig, axes = plt.subplots(4,4, figsize=(12,12))\nfor i, ax in enumerate(axes.reshape(-1)):\n    ax.imshow(img[:,:,1 + i])\nplt.show()\n\n# Check class imbalance\nprint(\"Aneurysm Present: 1 =\", train_df['Aneurysm Present'].mean()*100, \"%\")\n# Check modality distribution\nprint(train_df['Modality'].value_counts())\n# Check location-wise prevalence (critical for multi-label)\nlocations = [col for col in train_df.columns if 'Artery' in col or 'Communicating' in col]\nprint(train_df[locations].sum() / len(train_df))\n\ntrain_df.head()\n\ntrain_df['PatientAge'].describe()\n\ntrain_df['PatientAge'] = train_df['PatientAge'].astype(int)\nplt.figure(figsize=(10,6))\nsns.histplot(train_df['PatientAge'], bins=20, kde=False, color=sns.color_palette(\"rocket\")[4])  \nplt.xlabel('Patient Age')\nplt.ylabel('Count')\nplt.title('Distribution of Patient Age')\nplt.show()\n\nimport seaborn as sns\nimport matplotlib.pyplot as plt\n\n\n# Create cross-tabulation with proportions or counts\nctab = pd.crosstab(train_df['PatientSex'], train_df['Aneurysm Present'])\n\n\n# Plot grouped bar chart\nctab.plot(kind='bar', \n          color=sns.color_palette(\"pastel\"), \n          figsize=(8, 6), \n          width=0.8)\n\n\n# Labels and title\nplt.xlabel('Patient Sex')\nplt.ylabel('Count')\nplt.title('Aneurysm Presence by Patient Sex')\nplt.legend(title='Aneurysm Present', labels=['No', 'Yes'])\nplt.xticks(rotation=0)\n\n\n# Add value labels on bars (optional, improves readability)\nfor container in plt.gca().containers:\n    plt.bar_label(container, fmt='%d', padding=3)\n\n\nplt.tight_layout()\nplt.show()\n\n# Plot pie chart with counts shown on each slice\ntrain_df['Modality'].value_counts().plot(kind='pie', autopct='%d')\n\n\n# Optional: Improve layout and title\nplt.title('Distribution of Modality')\nplt.ylabel('')  # Hide the y-label (default is 'Modality' from pandas)\nplt.show()\n\nlabel_df.head()\n\nID_COL = 'SeriesInstanceUID'\n\n\nLABEL_COLS = [\n    'Left Infraclinoid Internal Carotid Artery',\n    'Right Infraclinoid Internal Carotid Artery',\n    'Left Supraclinoid Internal Carotid Artery',\n    'Right Supraclinoid Internal Carotid Artery',\n    'Left Middle Cerebral Artery',\n    'Right Middle Cerebral Artery',\n    'Anterior Communicating Artery',\n    'Left Anterior Cerebral Artery',\n    'Right Anterior Cerebral Artery',\n    'Left Posterior Communicating Artery',\n    'Right Posterior Communicating Artery',\n    'Basilar Tip',\n    'Other Posterior Circulation',\n    'Aneurysm Present',\n]\n\n\nDICOM_TAG_ALLOWLIST = [\n    'BitsAllocated', 'BitsStored', 'Columns', 'FrameOfReferenceUID', 'HighBit',\n    'ImageOrientationPatient', 'ImagePositionPatient', 'InstanceNumber', 'Modality',\n    'PatientID', 'PhotometricInterpretation', 'PixelRepresentation', 'PixelSpacing',\n    'PlanarConfiguration', 'RescaleIntercept', 'RescaleSlope', 'RescaleType', 'Rows',\n    'SOPClassUID', 'SOPInstanceUID', 'SamplesPerPixel', 'SliceThickness',\n    'SpacingBetweenSlices', 'StudyInstanceUID', 'TransferSyntaxUID',\n]\n\n\n# Model configuration\nTARGET_SIZE = (64, 64, 64)  # Reduced size for memory efficiency\nDEVICE = torch.device('cuda' if torch.cuda.is_available() else 'cpu')\n\n\nclass DICOMProcessor:\n    \"\"\"Process DICOM series into normalized 3D volumes\"\"\"\n    \n    def __init__(self, target_size: Tuple[int, int, int] = TARGET_SIZE):\n        self.target_size = target_size\n        self.scaler = StandardScaler()\n    \n    def load_dicom_series(self, series_path: str) -> np.ndarray:\n        \"\"\"Load and process a DICOM series into a 3D volume\"\"\"\n        try:\n            # Get all DICOM files\n            dicom_files = []\n            for root, _, files in os.walk(series_path):\n                for file in files:\n                    if file.endswith('.dcm'):\n                        dicom_files.append(os.path.join(root, file))\n            \n            if not dicom_files:\n                raise ValueError(f\"No DICOM files found in {series_path}\")\n            \n            # Load DICOMs \n            dicoms = []\n            # Initialize dicoms_with_position here (critical fix)\n            dicoms_with_position = []\n            for filepath in dicom_files:\n                try:\n                    ds = pydicom.dcmread(filepath, force=True)\n                    if hasattr(ds, 'PixelData'):\n                        # Extract Z-position (ImagePositionPatient[2]) with robust error handling\n                        z_pos = None\n                        if hasattr(ds, 'ImagePositionPatient') and len(ds.ImagePositionPatient) >= 3:\n                            try:\n                                z_pos = float(ds.ImagePositionPatient[2])\n                            except (ValueError, TypeError):\n                                # Handle case where ImagePositionPatient has non-numeric values\n                                pass\n                        dicoms_with_position.append((ds, filepath, z_pos))\n                except Exception as e:\n                    print(f\"Error reading {filepath}: {e}\")\n                    continue\n            \n            if not dicoms:\n                raise ValueError(f\"No valid DICOM files with pixel data in {series_path}\")\n            \n            # Sort by instance number\n            # dicoms.sort(key=lambda x: getattr(x[0], 'InstanceNumber', 0))\n            # CRITICAL FIX: Sort by physical Z-position, not InstanceNumber\n            # First, filter out any slices without valid Z-position\n            valid_slices = [(ds, filepath, z_pos) for ds, filepath, z_pos in dicoms_with_position if z_pos is not None]\n            \n            if not valid_slices:\n                # If no slices have valid Z-position, fall back to InstanceNumber as last resort\n                print(f\"Warning: No valid ImagePositionPatient found in {series_path}, using InstanceNumber\")\n                dicoms_with_position.sort(key=lambda x: getattr(x[0], 'InstanceNumber', 0))\n                dicoms = [(ds, filepath) for ds, filepath, _ in dicoms_with_position]\n            else:\n                # Primary method: Sort by physical Z-position\n                valid_slices.sort(key=lambda x: x[2])  # Sort by z_pos\n                dicoms = [(ds, filepath) for ds, filepath, _ in valid_slices]\n            \n            # Extract volume\n            volume_slices = []\n            for ds, _ in dicoms:\n                try:\n                    # Get pixel array\n                    pixel_array = ds.pixel_array.astype(np.float32)\n                    \n                    # Apply rescale if available\n                    if hasattr(ds, 'RescaleSlope') and hasattr(ds, 'RescaleIntercept'):\n                        slope = float(ds.RescaleSlope)\n                        intercept = float(ds.RescaleIntercept)\n                        pixel_array = pixel_array * slope + intercept\n                    \n                    volume_slices.append(pixel_array)\n                except Exception as e:\n                    print(f\"Error processing slice: {e}\")\n                    continue\n            \n            if not volume_slices:\n                raise ValueError(\"No valid slices extracted\")\n            \n            # Stack into 3D volume\n            volume = np.stack(volume_slices, axis=0)  # Shape: (depth, height, width)\n            \n            # Normalize and resize\n            volume = self.preprocess_volume(volume)\n            \n            return volume\n            \n        except Exception as e:\n            print(f\"Error processing series {series_path}: {e}\")\n            # Return zeros if processing fails\n            return np.zeros(self.target_size, dtype=np.float32)\n    \n    def preprocess_volume(self, volume: np.ndarray) -> np.ndarray:\n        \"\"\"Preprocess 3D volume: normalize, clip, resize\"\"\"\n        # Handle potential issues\n        if volume.size == 0:\n            return np.zeros(self.target_size, dtype=np.float32)\n        \n        # Clip extreme values (robust to outliers)\n        p1, p99 = np.percentile(volume, [1, 99])\n        volume = np.clip(volume, p1, p99)\n        \n        # Normalize to [0, 1]\n        volume_min, volume_max = volume.min(), volume.max()\n        if volume_max > volume_min:\n            volume = (volume - volume_min) / (volume_max - volume_min)\n        \n        # Resize to target size\n        if volume.shape != self.target_size:\n            zoom_factors = [\n                self.target_size[i] / volume.shape[i] for i in range(3)\n            ]\n            volume = ndimage.zoom(volume, zoom_factors, order=1)\n        \n        return volume.astype(np.float32)\n\n\nclass Simple3DCNN(nn.Module):\n    \"\"\"Lightweight 3D CNN for aneurysm detection\"\"\"\n    \n    def __init__(self, num_classes: int = len(LABEL_COLS)):\n        super(Simple3DCNN, self).__init__()\n        \n        # 3D Convolutional layers\n        self.conv1 = nn.Conv3d(1, 16, kernel_size=3, padding=1)\n        self.pool1 = nn.MaxPool3d(2)\n        self.conv2 = nn.Conv3d(16, 32, kernel_size=3, padding=1)\n        self.pool2 = nn.MaxPool3d(2)\n        self.conv3 = nn.Conv3d(32, 64, kernel_size=3, padding=1)\n        self.pool3 = nn.MaxPool3d(2)\n        self.conv4 = nn.Conv3d(64, 128, kernel_size=3, padding=1)\n        self.pool4 = nn.MaxPool3d(2)\n        \n        # Adaptive pooling to handle variable sizes\n        self.adaptive_pool = nn.AdaptiveAvgPool3d((2, 2, 2))\n        \n        # Fully connected layers\n        self.fc1 = nn.Linear(128 * 2 * 2 * 2, 256)\n        self.dropout1 = nn.Dropout(0.5)\n        self.fc2 = nn.Linear(256, 128)\n        self.dropout2 = nn.Dropout(0.3)\n        self.fc3 = nn.Linear(128, num_classes)\n        \n        # Batch normalization\n        self.bn1 = nn.BatchNorm3d(16)\n        self.bn2 = nn.BatchNorm3d(32)\n        self.bn3 = nn.BatchNorm3d(64)\n        self.bn4 = nn.BatchNorm3d(128)\n        \n    def forward(self, x):\n        # Input shape: (batch_size, 1, depth, height, width)\n        x = self.pool1(F.relu(self.bn1(self.conv1(x))))\n        x = self.pool2(F.relu(self.bn2(self.conv2(x))))\n        x = self.pool3(F.relu(self.bn3(self.conv3(x))))\n        x = self.pool4(F.relu(self.bn4(self.conv4(x))))\n        \n        # Adaptive pooling\n        x = self.adaptive_pool(x)\n        \n        # Flatten\n        x = x.view(x.size(0), -1)\n        \n        # Fully connected layers\n        x = F.relu(self.fc1(x))\n        x = self.dropout1(x)\n        x = F.relu(self.fc2(x))\n        x = self.dropout2(x)\n        x = self.fc3(x)\n        \n        return torch.sigmoid(x)\n\n\nclass AneurysmDataset(Dataset):\n    \"\"\"Dataset for loading training data\"\"\"\n    \n    def __init__(self, data_df: pd.DataFrame, series_dir: str, processor: DICOMProcessor):\n        self.data_df = data_df\n        self.series_dir = series_dir\n        self.processor = processor\n        \n    def __len__(self):\n        return len(self.data_df)\n    \n    def __getitem__(self, idx):\n        row = self.data_df.iloc[idx]\n        series_id = row[ID_COL]\n        \n        # Load volume\n        series_path = os.path.join(self.series_dir, series_id)\n        volume = self.processor.load_dicom_series(series_path)\n        \n        # Get labels\n        labels = row[LABEL_COLS].values.astype(np.float32)\n        \n        # Convert to tensor and add channel dimension\n        volume_tensor = torch.from_numpy(volume).unsqueeze(0)  # Add channel dim\n        labels_tensor = torch.from_numpy(labels)\n        \n        return volume_tensor, labels_tensor\n\n# Global model and processor\nmodel = None\nprocessor = None\n\n\ndef initialize_model():\n    \"\"\"Initialize model and processor (called once)\"\"\"\n    global model, processor\n    \n    if model is not None:\n        return\n    \n    print(\"Initializing model...\")\n    processor = DICOMProcessor(TARGET_SIZE)\n    model = Simple3DCNN(num_classes=len(LABEL_COLS))\n    \n    # Load pre-trained weights if available\n    try:\n        if os.path.exists('/kaggle/input/model_weights.pth'):\n            model.load_state_dict(torch.load('/kaggle/input/model_weights.pth', map_location='cpu'))\n            print(\"Loaded pre-trained weights\")\n        else:\n            print(\"No pre-trained weights found, using random initialization\")\n    except Exception as e:\n        print(f\"Error loading weights: {e}\")\n    \n    model.to(DEVICE)\n    model.eval()\n    print(f\"Model initialized on {DEVICE}\")\n\ndef train_model(series_dir: str = \"/kaggle/input/rsna-intracranial-aneurysm-detection/series\", num_epochs: int = 100, batch_size: int = 32):\n    \"\"\"Training function (for reference - would be run separately)\"\"\"\n    \n    # Load training data\n    train_df = pd.read_csv(\"/kaggle/input/rsna-intracranial-aneurysm-detection/train.csv\")\n    \n    # Initialize components\n    processor = DICOMProcessor(TARGET_SIZE)\n    model = Simple3DCNN(num_classes=len(LABEL_COLS))\n    model.to(DEVICE)\n    \n    # Create dataset and dataloader\n    dataset = AneurysmDataset(train_df, series_dir, processor)\n    dataloader = DataLoader(dataset, batch_size=batch_size, shuffle=True, num_workers=0)\n    \n    # Loss and optimizer\n    criterion = nn.BCELoss()\n    optimizer = optim.Adam(model.parameters(), lr=0.001, weight_decay=1e-5)\n    scheduler = optim.lr_scheduler.ReduceLROnPlateau(optimizer, mode='min', patience=2)\n    \n    # Training loop\n    model.train()\n    for epoch in range(num_epochs):\n        total_loss = 0\n        num_batches = 0\n        \n        for batch_idx, (volumes, labels) in enumerate(dataloader):\n            volumes, labels = volumes.to(DEVICE), labels.to(DEVICE)\n            \n            optimizer.zero_grad()\n            outputs = model(volumes)\n            loss = criterion(outputs, labels)\n            loss.backward()\n            optimizer.step()\n            \n            total_loss += loss.item()\n            num_batches += 1\n            \n            if batch_idx % 10 == 0:\n                print(f'Epoch {epoch+1}/{num_epochs}, Batch {batch_idx+1}, Loss: {loss.item():.4f}')\n        \n        avg_loss = total_loss / num_batches\n        scheduler.step(avg_loss)\n        print(f'Epoch {epoch+1}/{num_epochs} completed, Average Loss: {avg_loss:.4f}')\n    \n    # Save model\n    torch.save(model.state_dict(), 'model_weights.pth')\n    print(\"Model saved to model_weights.pth\")\n\n# train_model()\n\ndef predict(series_path: str) -> pl.DataFrame:\n    \"\"\"Make prediction for a single series\"\"\"\n    \n    # Initialize model on first call\n    initialize_model()\n    \n    series_id = os.path.basename(series_path)\n    \n    try:\n        # Process the DICOM series\n        volume = processor.load_dicom_series(series_path)\n        \n        # Convert to tensor and add batch dimension\n        volume_tensor = torch.from_numpy(volume).unsqueeze(0).unsqueeze(0)  # (1, 1, D, H, W)\n        volume_tensor = volume_tensor.to(DEVICE)\n        \n        # Make prediction\n        with torch.no_grad():\n            predictions = model(volume_tensor)\n            predictions = predictions.cpu().numpy().flatten()\n        \n        # Create result DataFrame\n        result_data = [[series_id] + predictions.tolist()]\n        result_df = pl.DataFrame(\n            data=result_data,\n            schema=[ID_COL] + LABEL_COLS,\n            orient='row'\n        )\n        \n        # Clean up memory\n        del volume_tensor\n        torch.cuda.empty_cache()\n        gc.collect()\n        \n    except Exception as e:\n        print(f\"Error predicting for series {series_id}: {e}\")\n        # Return baseline predictions (0.5 for all classes)\n        result_data = [[series_id] + [0.5] * len(LABEL_COLS)]\n        result_df = pl.DataFrame(\n            data=result_data,\n            schema=[ID_COL] + LABEL_COLS,\n            orient='row'\n        )\n    \n    # Mandatory cleanup\n    shutil.rmtree('/kaggle/shared', ignore_errors=True)\n    \n    return result_df.drop(ID_COL)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-28T23:49:01.95115Z","iopub.execute_input":"2025-09-28T23:49:01.951369Z","iopub.status.idle":"2025-09-28T23:49:13.957331Z","shell.execute_reply.started":"2025-09-28T23:49:01.951346Z","shell.execute_reply":"2025-09-28T23:49:13.956518Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"inference_server = kaggle_evaluation.rsna_inference_server.RSNAInferenceServer(predict)\n\nif os.getenv('KAGGLE_IS_COMPETITION_RERUN'):\n    inference_server.serve()\nelse:\n    inference_server.run_local_gateway()\n    display(pl.read_parquet('/kaggle/working/submission.parquet'))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-09-28T23:49:13.958215Z","iopub.execute_input":"2025-09-28T23:49:13.958505Z","iopub.status.idle":"2025-09-28T23:49:28.862405Z","shell.execute_reply.started":"2025-09-28T23:49:13.958479Z","shell.execute_reply":"2025-09-28T23:49:28.86177Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}