Dice similarity index between two 3D array from rtstructure dicom file

Question

hello all i created a script that evaluates dice similarity index (DSC) but i'm not sure of the script. after selecting structures i want to calculate their DSC i define two funtion the first one to extract the array: \python\ code:

def extract_coordinates(rtstruct, target_name):
    for item in rtstruct.StructureSetROISequence:
        roi_name = item.ROIName.strip()
        roi_number = item.ROINumber
        if roi_name.lower() == target_name.lower():
            
            referenced_roi_number = rtstruct.ROIContourSequence[roi_number -1].ReferencedROINumber
            coordinate_array = []
            for dataset in rtstruct.ROIContourSequence[roi_number -1].ContourSequence:
                coordinate_triplets = dataset.ContourData
                for i in range(0, len(coordinate_triplets), 3):
                    x, y, z = map(float, coordinate_triplets[i:i + 3])
                    coordinate_array.append((x, y, z))
            return coordinate_array
    return None

coordinate_array1 = extract_coordinates(rtstruct1, target_name1)
coordinate_array2 = extract_coordinates(rtstruct2, target_name2)

the lenth of the list is about 32 objects. if i select two final two elements i got an error on line referenced_roi_number = rtstruct.ROIContourSequence[roi_number -1].ReferencedROINumber-> out of index, moreover changing [roi_number-1] to [roi_number-10] the DSC will change but i0m able to use the final two item. the second function is about how i calculate the DSC, but if i swap the PTV (contour1, while ctv contour2) with CTV (that now became contour1 with PTV contour2), the DSC will change from 0.2853 to 1.7147 (2-1.7147= 0.2853) but i need always the same result not depending on the order of structures.

if coordinate_array1 is not None and coordinate_array2 is not None:
    def calculate_dice_similarity(coordinate_array1, coordinate_array2):
        contour1 = set(coordinate_array1)  # Converti la lista in un insieme
        contour2 = set(coordinate_array2)  # Converti la lista in un insieme
    
        intersection = len(np.logical_and(contour1, contour2))
        union = len(np.union1d(contour1, contour2))
    
        dice_similarity = (2.0 * intersection) / (len(contour1) + len(contour2))
        return dice_similarity
        

    # Evaluation of  Dice Similarity Index
    dice_similarity = calculate_dice_similarity(coordinate_array1, coordinate_array2)
    print(f"Dice Similarity Index between '{target_name1}' e '{target_name2}': {dice_similarity:.4f}")
else:
    print(f"Struture '{target_name1}' o '{target_name2}' not found

i've got no idea for resolving them. i tried to access the list of roicontoursequence and i find out that is a list of 32 elements but i don't know ho to fix the problem. The same for DSC, the index does not depend on order of factor. thank for all

i've used the suggetion of the answer by making this correction to my code, here i reported both the modified extraction function and the DSC calculation function:

ef extract_coordinates(rtstruct, target_name):
for item in rtstruct.StructureSetROISequence:
    roi_to_name = {roi.ROINumber: roi.ROIName for roi in rtstruct.StructureSetROISequence}
    # get contour index->ROIName
    ctr_index_to_name = {i: roi_to_name[ctr.ReferencedROINumber] for i, ctr in enumerate(rtstruct.ROIContourSequence)}
    # Reverse to get name->index
    name_to_ctr_index = {v: k for k,v in ctr_index_to_name.items()}
    referenced_roi_number = rtstruct.ROIContourSequence[name_to_ctr_index[target_name]].ReferencedROINumber
    coordinate_array = []
    
        
        
    for dataset in rtstruct.ROIContourSequence[name_to_ctr_index[target_name]].ContourSequence:
            coordinate_triplets = dataset.ContourData
            for i in range(0, len(coordinate_triplets), 3):
                x, y, z = map(float, coordinate_triplets[i:i + 3])
                coordinate_array.append((x, y, z))
                return coordinate_array
return None

the DSC calculation function remained:

coordinate_array2 = extract_coordinates(rtstruct2, target_name2)

if coordinate_array1 is not None and coordinate_array2 is not None:
    def calculate_dice_similarity(coordinate_array1, coordinate_array2):
        contour1 = set(coordinate_array1)  # Converti la lista in un insieme
        contour2 = set(coordinate_array2)  # Converti la lista in un insieme
    
        intersection = len(np.logical_and(contour1, contour2))
        union = len(np.union1d(contour1, contour2))
    
        dice_similarity = (2.0 * intersection) / (len(contour1) + len(contour2))
        return dice_similarity
        

    # Calcola il Dice Similarity Index
    dice_similarity = calculate_dice_similarity(coordinate_array1, coordinate_array2)
    print(f"Dice Similarity Index tra '{target_name1}' e '{target_name2}': {dice_similarity:.4f}")
else:
    print(f"Structure '{target_name1}' o '{target_name2}' not found nel file RTSTRUCT.")

but this time the output is always 1.0 since i followed the last comment suggestion i tried to define internal polygon points redifining the calculation with this method:

coordinate_array1 = extract_coordinates(rtstruct1, target_name1)

coordinate_array2 = extract_coordinates(rtstruct2, target_name2)

if coordinate_array1 is not None and coordinate_array2 is not None: def calculate_dice_similarity(coordinate_array1, coordinate_array2): # Convert the coordinates to Shapely polygons polygon1 = Polygon(coordinate_array1) polygon2 = Polygon(coordinate_array2)

# Generate a set of all points within the bounding boxes of both polygons
 bounding_box1 = polygon1.bounds
 bounding_box2 = polygon2.bounds
 min_x = min(bounding_box1[0], bounding_box2[0])
 max_x = max(bounding_box1[2], bounding_box2[2])
 min_y = min(bounding_box1[1], bounding_box2[1])
 max_y = max(bounding_box1[3], bounding_box2[3])

 points_inside = set()
 for x in range(int(min_x), int(max_x) + 1):
     for y in range(int(min_y), int(max_y) + 1):
         point = Point(x, y)
         if polygon1.contains(point) or polygon2.contains(point):
             points_inside.add(point)

# Calculate the Dice similarity based on the number of points inside
 intersection = len(points_inside.intersection(polygon1).intersection(polygon2))
 union = len(points_inside.intersection(polygon1).union(polygon2))

 dice_similarity = (2.0 * intersection) / (len(polygon1) + len(polygon2))
 return dice_similarity

 print(f"Dice Similarity Index tra '{target_name1}' e '{target_name2}': {dice_similarity:.4f}")

else: print(f"Struttura '{target_name1}' o '{target_name2}' non trovata nel file RTSTRUCT.")

DSC= calculate_dice_similarity(coordinate_array1, coordinate_array2) but i have this error message: ` File D:\spyder_env\Lib\site-packages\spyder_kernels\py3compat.py:356 in compat_exec exec(code, globals, locals)

File d:\scipts_homemade\dice_test_2.py:128 DSC= calculate_dice_similarity(coordinate_array1, coordinate_array2)

File d:\scipts_homemade\dice_test_2.py:99 in calculate_dice_similarity polygon1 = Polygon(coordinate_array1)

File D:\spyder_env\Lib\site-packages\shapely\geometry\polygon.py:230 in new shell = LinearRing(shell)

File D:\spyder_env\Lib\site-packages\shapely\geometry\polygon.py:104 in new geom = shapely.linearrings(coordinates)

File D:\spyder_env\Lib\site-packages\shapely\decorators.py:77 in wrapped return func(*args, **kwargs)

File D:\spyder_env\Lib\site-packages\shapely\creation.py:171 in linearrings return lib.linearrings(coords, out=out, **kwargs)

ValueError: A linearring requires at least 4 coordinates.

#so i printed the array, for example the one sample print(coordinate_array1) [(-12.88, -211.57, -339.0)]` how can i resolve this?

Answer 1

Given that you are getting an IndexError, without further information I am guessing that the ROI number does not represent the correct index into the contour array. They are often stored that way but it is not guaranteed.

I would suggest something like the following to get a mapping of target_name to ROIContour index:

roi_to_name = {roi.ROINumber: roi.ROIName for roi in ds.StructureSetROISequence}
# get ROIname -> contour index
name_to_ctr_index = {roi_to_name[ctr.ReferencedROINumber]: i for i, ctr in enumerate(ds.ROIContourSequence)}

Then you can jump straight to:

for dataset in rtstruct.ROIContourSequence[name_to_ctr_index[target_name]].ContourSequence:
   ...