Skip to content

Commit

Permalink
added agatson option to aortic calcium
Browse files Browse the repository at this point in the history
  • Loading branch information
malteekj committed Jan 31, 2024
1 parent 57170db commit de462a8
Show file tree
Hide file tree
Showing 4 changed files with 127 additions and 21 deletions.
12 changes: 9 additions & 3 deletions bin/C2C
Original file line number Diff line number Diff line change
Expand Up @@ -128,7 +128,8 @@ def AorticCalciumPipelineBuilder(path, args):
aortic_calcium.AorticCalciumMetrics(),
aortic_calcium_visualization.AorticCalciumVisualizer(),
aortic_calcium_visualization.AorticCalciumPrinter(),
]
],
args=args
)
return pipeline

Expand Down Expand Up @@ -201,8 +202,13 @@ def argument_parser():
"liver_spleen_pancreas", parents=[base_parser]
)

aortic_calcium = subparsers.add_parser("aortic_calcium", parents=[base_parser])

aortic_calcium = subparsers.add_parser(
"aortic_calcium", parents=[base_parser])

aortic_calcium.add_argument(
"--threshold", default="adaptive", type=str
)

contrast_phase_parser = subparsers.add_parser(
"contrast_phase", parents=[base_parser]
)
Expand Down
79 changes: 64 additions & 15 deletions comp2comp/aortic_calcium/aortic_calcium.py
Original file line number Diff line number Diff line change
Expand Up @@ -126,7 +126,8 @@ def __call__(self, inference_pipeline):
spine_mask = inference_pipeline.spine_segmentation.get_fdata() > 0

calcification_results = self.detectCalcifications(
ct, aorta_mask, exclude_mask=spine_mask, remove_size=3, return_dilated_mask=True
ct, aorta_mask, exclude_mask=spine_mask, remove_size=3, return_dilated_mask=True,
threshold=inference_pipeline.args.threshold
)

inference_pipeline.calc_mask = calcification_results['calc_mask']
Expand Down Expand Up @@ -182,6 +183,7 @@ def detectCalcifications(
exclude_center_aorta=True,
return_eroded_aorta=False,
aorta_erode_iteration=6,
threshold = 'adaptive',
):
"""
Function that takes in a CT image and aorta segmentation (and optionally volumes to use
Expand Down Expand Up @@ -222,6 +224,9 @@ def detectCalcifications(
Return the eroded center aorta. Defaults to False.
aorta_erode_iteration (int, optional):
Number of iterations for the strcturing element. Defaults to 6.
threshold: (str, int):
Can either be 'adaptive', 'agatson', or int. Choosing 'agatson'
Will mean a threshold of 130 HU.
Returns:
results: array of only the mask is returned, or dict if other volumes are also returned.
Expand Down Expand Up @@ -268,16 +273,25 @@ def slicedDilationOrErosion(input_mask, struct, num_iteration, operation):

biggest_idx = np.argmax(aorta_vols) + 1
aorta_mask[labelled_aorta != biggest_idx] = 0

# Get aortic CT point to set adaptive threshold
aorta_ct_points = ct[aorta_mask == 1]

# equal to one standard deviation to the left of the curve
quant = 0.158
quantile_median_dist = np.median(aorta_ct_points) - np.quantile(
aorta_ct_points, q=quant
)
calc_thres = np.median(aorta_ct_points) + quantile_median_dist * num_std

### Choose the threshold ###
if threshold == 'adaptive':
# Get aortic CT point to set adaptive threshold
aorta_ct_points = ct[aorta_mask == 1]

# equal to one standard deviation to the left of the curve
quant = 0.158
quantile_median_dist = np.median(aorta_ct_points) - np.quantile(
aorta_ct_points, q=quant
)
calc_thres = np.median(aorta_ct_points) + quantile_median_dist * num_std
elif threshold == 'agatson':
calc_thres = 130
elif isinstance(threshold, int):
calc_thres = threshold
else:
raise ValueError('Error in threshold value for aortic calcium segmentaiton. \
Should be \'adaptive\', \'agatson\' or int, but got: ' + str(threshold))

t0 = time.time()

Expand Down Expand Up @@ -448,15 +462,17 @@ def __call__(self, inference_pipeline):

ct_full = inference_pipeline.medical_volume.get_fdata()

for i in range(2):
for i in range(len(region_names)):
# count statistics for individual calcifications
if i == 0:
labelled_calc, num_lesions = ndimage.label(calc_mask[:,:,:sep_plane])
calc_mask_region = calc_mask[:,:,:sep_plane]
ct = ct_full[:,:,:sep_plane]
elif i == 1:
labelled_calc, num_lesions = ndimage.label(calc_mask[:,:,sep_plane:])
calc_mask_region = calc_mask[:,:,sep_plane:]
ct = ct_full[:,:,sep_plane:]

labelled_calc, num_lesions = ndimage.label(calc_mask_region)

metrics = {
"volume": [],
"mean_hu": [],
Expand All @@ -481,9 +497,42 @@ def __call__(self, inference_pipeline):
metrics["volume_total"] = calc_vol

metrics["num_calc"] = len(metrics["volume"])


if inference_pipeline.args.threshold == 'agatson':
metrics["agatson_score"] = self.CalculateAgatsonScore(calc_mask_region, ct, inference_pipeline.pix_dims)

all_regions[region_names[i]] = metrics

inference_pipeline.metrics = all_regions

return {}

def CalculateAgatsonScore(self, calc_mask_region, ct, pix_dims):
def get_hu_factor(max_hu):
# if max_hu ><
if max_hu < 200:
factor = 1
elif 200 <= max_hu < 300:
factor = 2
elif 300 <= max_hu < 400:
factor = 3
elif max_hu > 400:
factor = 4

return factor

# dims are in mm here
area_per_pixel = pix_dims[0] * pix_dims[1]
agatson = 0

for i in range(calc_mask_region.shape[2]):
tmp_slice = calc_mask_region[:,:,i]
tmp_ct_slice = ct[:,:,i]

labelled_calc, num_lesions = ndimage.label(tmp_slice)

for j in range(1, num_lesions + 1):
tmp_mask = labelled_calc == j
agatson += tmp_mask.sum() * area_per_pixel * get_hu_factor(tmp_ct_slice[tmp_mask].max())

return agatson
54 changes: 52 additions & 2 deletions comp2comp/aortic_calcium/aortic_calcium_visualization.py
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,7 @@ def __call__(self, inference_pipeline):
):
f.write("{},{:.1f},{:.1f},{:.1f}\n".format(vol, mean, median, max))


# Write total results
with open(
os.path.join(
Expand All @@ -68,8 +69,49 @@ def __call__(self, inference_pipeline):
f.write(region + ',\n')

f.write("Total number,{}\n".format(metrics["num_calc"]))
f.write("Total volume (cm^3),{}\n".format(metrics["volume_total"]))
f.write("Threshold (HU),{}\n".format(inference_pipeline.calcium_threshold))
f.write("Total volume (cm^3),{:.1f}\n".format(metrics["volume_total"]))
f.write("Threshold (HU),{:.1f}\n".format(inference_pipeline.calcium_threshold))

f.write("{},{:.1f}+/-{:.1f}\n".format(
"Mean HU",
np.mean(metrics["mean_hu"]),
np.std(metrics["mean_hu"]),
)
)
f.write("{},{:.1f}+/-{:.1f}\n".format(
"Median HU",
np.mean(metrics["median_hu"]),
np.std(metrics["median_hu"]),
)
)
f.write("{},{:.1f}+/-{:.1f}\n".format(
"Max HU",
np.mean(metrics["max_hu"]),
np.std(metrics["max_hu"]),
)
)
f.write("{},{:.1f}+/-{:.1f}\n".format(
"Mean volume (cm³):",
np.mean(metrics["volume"]),
np.std(metrics["volume"]),
)
)
f.write("{},{:.1f}\n".format(
"Median volume (cm³)", np.median(metrics["volume"])
)
)
f.write("{},{:.1f}\n".format(
"Max volume (cm³)", np.max(metrics["volume"])
)
)
f.write("{},{:.1f}\n".format(
"Min volume (cm³):", np.min(metrics["volume"])
)
)

if inference_pipeline.args.threshold == 'agatson':
f.write("Agatson score,{:.1f}\n".format(metrics["agatson_score"]))


distance = 25
print("\n")
Expand Down Expand Up @@ -139,7 +181,15 @@ def __call__(self, inference_pipeline):
"Threshold (HU):", distance, inference_pipeline.calcium_threshold
)
)
if inference_pipeline.args.threshold == 'agatson':
print(
"{:<{}}{:.1f}".format(
"Agatson score:", distance, metrics["agatson_score"]
)
)


print("\n")

return {}

3 changes: 2 additions & 1 deletion comp2comp/inference_pipeline.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,8 +14,9 @@
class InferencePipeline(InferenceClass):
"""Inference pipeline."""

def __init__(self, inference_classes: List = None, config: Dict = None):
def __init__(self, inference_classes: List = None, config: Dict = None, args: Dict = None):
self.config = config
self.args = args
# assign values from config to attributes
if self.config is not None:
for key, value in self.config.items():
Expand Down

0 comments on commit de462a8

Please sign in to comment.