# Dice Coefficient Calculation in Python

**URL:** <https://discourse.itk.org/t/dice-coefficient-calculation-in-python/7628>\
**Category:** Beginner Questions\
**Tags:** github, dicom, simpleitk, python\
**Created:** [September 5, 2025, 11:35am UTC](https://discourse.itk.org/t/dice-coefficient-calculation-in-python/7628 "2025-09-05T11:35:47Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![blue\_sky](https://discourse.itk.org/letter_avatar_proxy/v4/letter/b/dbc845/32.png) [@blue\_sky](https://discourse.itk.org/u/blue_sky)\
**Post date:** [September 5, 2025, 11:35am UTC](https://discourse.itk.org/t/dice-coefficient-calculation-in-python/7628/1 "2025-09-05T11:35:47Z")

</div>

Hello everyone,

I’m trying to compute the Dice coefficient using a Python script.

To implement this, I have used ideas and code from the following links:

“[SimpleITK-Notebooks/Python/34\_Segmentation\_Evaluation.ipynb at main · InsightSoftwareConsortium/SimpleITK-Notebooks · GitHub](https://github.com/InsightSoftwareConsortium/SimpleITK-Notebooks/blob/main/Python/34_Segmentation_Evaluation.ipynb)”

“[SlicerRT/SegmentComparison at master · SlicerRt/SlicerRT · GitHub](https://github.com/SlicerRt/SlicerRT/tree/master/SegmentComparison)”

However, as shown in the attached screenshot, there is a slight difference between the Dice coefficient values computed by my code and those reported by 3D Slicer  
(my results:  
Dice for Rectum=0.628, Bladder= 0.749)

Could you please help me understand why there is this discrepancy and how I can resolve it so that my results match Slicer’s output?

Thanks a lot

my code:

> import SimpleITK as sitk  
> import numpy as np  
> import pydicom  
> import logging  
> import os  
> from skimage.draw import polygon  
> from scipy.spatial.distance import directed\_hausdorff  
> from scipy.ndimage import binary\_dilation  
> logging.basicConfig(level=logging.INFO, format=‘%(asctime)s - %(levelname)s - %(message)s’)  
> def read\_dicom\_series(directory):  
> if not os.path.isdir(directory): raise FileNotFoundError(f"Directory not found: {directory}“)  
> series\_IDs = sitk.ImageSeriesReader.GetGDCMSeriesIDs(directory)  
> if not series\_IDs: raise ValueError(f"No DICOM series found in: {directory}”)  
> series\_files = sitk.ImageSeriesReader.GetGDCMSeriesFileNames(directory, series\_IDs[0])  
> reader = sitk.ImageSeriesReader(); reader.SetFileNames(series\_files)  
> return reader.Execute()  
> def list\_rois(rtstruct\_path):  
> ds = pydicom.dcmread(rtstruct\_path)  
> return [roi.ROIName for roi in ds.StructureSetROISequence] if hasattr(ds, ‘StructureSetROISequence’) else []  
> def extract\_mask\_from\_rtstruct(rtstruct\_path, roi\_name, reference\_ct, output\_nifti=None):  
> ds = pydicom.dcmread(rtstruct\_path)  
> roi\_number = next((i+1 for i, roi in enumerate(ds.StructureSetROISequence) if roi.ROIName.lower() == roi\_name.lower()), None)  
> if not roi\_number:  
> logging.warning(f"ROI ‘{roi\_name}’ not found in {rtstruct\_path}.“)  
> return None  
> contour\_seq = next((rc.ContourSequence for rc in ds.ROIContourSequence if rc.ReferencedROINumber == roi\_number), None)  
> if not contour\_seq:  
> logging.warning(f"No contours found for ‘{roi\_name}’ in {rtstruct\_path}.”)  
> return None  
> mask\_array = np.zeros(reference\_ct.GetSize()[::-1], dtype=np.uint8) # z, y, x  
> has\_contours = False  
> for contour in contour\_seq:  
> points = np.array(contour.ContourData).reshape(-1, 3)  
> if len(points) \< 3: continue  
> has\_contours = True  
> z\_pos = points[0, 2]  
> slice\_index = int(round(reference\_ct.TransformPhysicalPointToContinuousIndex((0, 0, z\_pos))[2]))  
> if not (0 \<= slice\_index \< mask\_array.shape[0]):  
> logging.warning(f"Slice index {slice\_index} out of bounds - skipping.“)  
> continue  
> rows, cols = [], []  
> for p in points:  
> cont\_idx = reference\_ct.TransformPhysicalPointToContinuousIndex(p)  
> rows.append(cont\_idx[1])  
> cols.append(cont\_idx[0])  
> rr, cc = polygon(rows, cols, shape=mask\_array.shape[1:])  
> mask\_array[slice\_index, rr, cc] = 1  
> if not has\_contours:  
> logging.warning(f"No valid contours for ‘{roi\_name}’ - returning None.”)  
> return None  
> ####  
> mask\_array = binary\_dilation(mask\_array, structure=np.ones((3,3,3)))  
> mask = sitk.GetImageFromArray(mask\_array.astype(np.uint8))  
> mask = sitk.BinaryFillhole(mask) # Fill holes like Slicer  
> mask.CopyInformation(reference\_ct)  
> if output\_nifti:  
> sitk.WriteImage(mask, output\_nifti)  
> logging.info(f"Mask saved to: {output\_nifti}“)  
> return mask  
> def compute\_dice(fixed\_mask, moving\_mask):  
> overlap\_filter = sitk.LabelOverlapMeasuresImageFilter()  
> overlap\_filter.Execute(fixed\_mask, moving\_mask)  
> return overlap\_filter.GetDiceCoefficient()  
> def compute\_hausdorff\_distance(fixed\_mask, moving\_mask):  
> fixed\_points = np.argwhere(sitk.GetArrayFromImage(fixed\_mask) \> 0)  
> moving\_points = np.argwhere(sitk.GetArrayFromImage(moving\_mask) \> 0)  
> if fixed\_points.size == 0 or moving\_points.size == 0:  
> return float(‘inf’)  
> dist1 = directed\_hausdorff(fixed\_points, moving\_points)[0]  
> dist2 = directed\_hausdorff(moving\_points, fixed\_points)[0]  
> return max(dist1, dist2)  
> def compute\_overlap\_fraction(fixed\_mask, moving\_mask):  
> intersection = np.logical\_and(sitk.GetArrayFromImage(fixed\_mask), sitk.GetArrayFromImage(moving\_mask)).sum()  
> return intersection / (sitk.GetArrayFromImage(fixed\_mask).sum() + 1e-6)  
> def main():  
> fixed\_ct\_dir = r’D:\fr4\CT\_Unified’  
> fixed\_rtstruct\_path = r’D:\fr4\STRUCT\_Unified\RTS\_fr4.dcm’  
> transformed\_rtstruct\_path = r’D:\fr2\STRUCT\_Unified\RTS\_fr2.dcm’  
> fixed\_ct = read\_dicom\_series(fixed\_ct\_dir)  
> fixed\_rois = list\_rois(fixed\_rtstruct\_path)  
> transformed\_rois = list\_rois(transformed\_rtstruct\_path)  
> logging.info(”\nAvailable ROIs in fixed RTSTRUCT:“)  
> for roi in fixed\_rois: logging.info(f” - {roi}“)  
> logging.info(”\nAvailable ROIs in transformed RTSTRUCT:“)  
> for roi in transformed\_rois: logging.info(f” - {roi}“)  
> organs = [‘Bladder’, ‘Rectum’]  
> for organ in organs:  
> try:  
> fixed\_roi = input(f"Enter ROI name for {organ} in fixed RTSTRUCT (from list above): “).strip()  
> transformed\_roi = input(f"Enter ROI name for {organ} in transformed RTSTRUCT (from list above): “).strip()  
> fixed\_mask = extract\_mask\_from\_rtstruct(fixed\_rtstruct\_path, fixed\_roi, fixed\_ct, f"fixed\_{organ}_mask.nii.gz")  
> transformed\_mask = extract\_mask\_from\_rtstruct(transformed\_rtstruct\_path, transformed\_roi, fixed\_ct, f"transformed_{organ}\_mask.nii.gz”)  
> if fixed\_mask is None or transformed\_mask is None:  
> logging.warning(f"Skipping metrics for ‘{organ}’ due to missing mask.”)  
> continue  
> resampler = sitk.ResampleImageFilter()  
> resampler.SetReferenceImage(fixed\_mask)  
> resampler.SetInterpolator(sitk.sitkNearestNeighbor)  
> transformed\_mask = resampler.Execute(transformed\_mask)  
> dice = compute\_dice(fixed\_mask, transformed\_mask)  
> hausdorff = compute\_hausdorff\_distance(fixed\_mask, transformed\_mask)  
> overlap = compute\_overlap\_fraction(fixed\_mask, transformed\_mask)  
> logging.info(f”\nMetrics for {organ} (fixed: ‘{fixed\_roi}’, transformed: ‘{transformed\_roi}’):“)  
> logging.info(f” Dice Coefficient: {dice:.4f}“)  
> logging.info(f” Hausdorff Distance: {hausdorff:.4f} mm")  
> logging.info(f" Overlap Fraction: {overlap:.4f}“)  
> except Exception as e:  
> logging.error(f"Error computing metrics for {organ}: {e}”)  
> if **name** == “ **main** ”:  
> main()

 ![Dice_1](https://discourse.itk.org/uploads/default/original/2X/1/137f410fac9eb857af4a8b14073fb98b33dfd1d5.png)

---

<div class="post-metadata">

**Author:** ![zivy](https://discourse.itk.org/user_avatar/discourse.itk.org/zivy/32/1726_2.png) [@zivy](https://discourse.itk.org/u/zivy)\
**Post date:** [September 5, 2025, 1:11pm UTC](https://discourse.itk.org/t/dice-coefficient-calculation-in-python/7628/2 "2025-09-05T13:11:36Z")

</div>

Hello @blue_sky,

I suspect this is due to differences in the way the RTSRUCT contour information is converted to masks.

In the provided code the conversion is done via skimage.draw.polygon (generate coordinates of pixels inside a polygon). I don’t know how it is done in Slicer but it is done differently.

To see if this hunch is valid, convert an RTSRUCT to a binary image using your code. Load this RTSTRUCT into Slicer, load the binary image you created in Slicer and compute the Dice coefficient between them.

You could do this the other way round too. Load RTSRUCT into Slicer, export it as binary image. Convert the RTSRUCT to a binary image using your code. Read both images using SimpleITK and compute the Dice coefficient.

As a side note, when sharing code, please ensure that the indentation is correct, makes the code more readable.

---

<div class="post-metadata">

**Author:** ![shai-ikko](https://discourse.itk.org/user_avatar/discourse.itk.org/shai-ikko/32/2842_2.png) [@shai-ikko](https://discourse.itk.org/u/shai-ikko)\
**Post date:** [September 11, 2025, 9:01am UTC](https://discourse.itk.org/t/dice-coefficient-calculation-in-python/7628/3 "2025-09-11T09:01:10Z")

</div>

> [@zivy](#):
>
> As a side note, when sharing code, please ensure that the indentation is correct, makes the code more readable.

@blue_sky use code blocks instead of quoting. If you place three backticks (the ` sign) you get a block that’s in monospace. On the right you can even choose the programming language for syntax highlighting.

```python
import os
import whatever

def f(x: int) -> float
    return x + 0.5 if "yes" else x + 0.1

```