Segment CT Volumes With Several Readers and STAPLE

Medical imaging researchers who segment CT scans in 3D Slicer or ITK-SNAP can follow this tutorial to give the same scans to three readers in Potato 2.13.1, read the readers' 3D Dice agreement in Potato, export each reader's masks as NIfTI label volumes, and build a consensus mask with STAPLE. STAPLE (Simultaneous Truth and Performance Level Estimation) estimates a consensus segmentation from several readers' masks, together with each reader's sensitivity and specificity against that consensus. The Dice coefficient measures the overlap of two masks, from 0 for no shared voxels to 1 for identical masks.

Potato 2.13.1 computes Dice between readers and, with three or more readers, STAPLE's sensitivity and specificity for each reader, and Potato does not export a consensus mask. The consensus mask in this tutorial comes from SimpleITK on the masks Potato exports. The worked example segments the spleen on three CT scans from the Medical Segmentation Decathlon. The numbers on this page come from a run on those scans with Potato 2.12.1, and the agreement report was run again on the saved masks with Potato 2.13.1. Potato is research software with no clinical or regulatory clearance, and nothing on this page is meant for patient care.

The scans

The Medical Segmentation Decathlon (MSD) released ten segmentation tasks for benchmarking (Antonelli et al., 2022). Task09 is the spleen, with 61 CT volumes from Memorial Sloan Kettering Cancer Center. The 41 training volumes come with an expert spleen mask, and the 20 test volumes come without one. The MSD site and the task's dataset.json both give the license as CC BY-SA 4.0. Under CC BY-SA 4.0, you can share and adapt the scans with attribution and must share adaptations under the same license.

Task09 downloads as one 1.6 GB archive. This tutorial uses three of the smallest training scans, about 30 MB together, and keeps the expert masks aside as an answer key for checking the consensus later:

curl -O https://msd-for-monai.s3-us-west-2.amazonaws.com/Task09_Spleen.tar
tar -xf Task09_Spleen.tar \
  Task09_Spleen/imagesTr/spleen_3.nii.gz Task09_Spleen/labelsTr/spleen_3.nii.gz \
  Task09_Spleen/imagesTr/spleen_8.nii.gz Task09_Spleen/labelsTr/spleen_8.nii.gz \
  Task09_Spleen/imagesTr/spleen_56.nii.gz Task09_Spleen/labelsTr/spleen_56.nii.gz

Each scan is 512 × 512 voxels per slice, with slices from 5.0 to 7.5 mm thick, and the spleen spans between 9 and 19 axial slices:

Installing Potato with volume support

The volume viewer reads NIfTI through nibabel and DICOM through pydicom, and writes DICOM SEG through highdicom. The volumes extra installs the three libraries:

pip install 'potato-annotation[volumes]==2.13.1'

Potato serves the viewer's scripts and styles from the installed package. A headless browser on the annotation page in the test run made no request to any other host. Once the packages and scans are on the machine, the whole tutorial runs without a network connection, and Air-Gapped Deployment covers that setup.

The project folder

Put the three scans in a media/ folder and list one case per line in a data file:

spleen-readers/
├── config.yaml
├── data/
│   └── cases.jsonl
└── media/
    ├── spleen_3.nii.gz
    ├── spleen_8.nii.gz
    └── spleen_56.nii.gz

volume names the scan relative to media/, and text is the caption shown above the viewer:

{"id": "spleen_3", "volume": "spleen_3.nii.gz", "text": "MSD Task09 spleen_3"}
{"id": "spleen_8", "volume": "spleen_8.nii.gz", "text": "MSD Task09 spleen_8"}
{"id": "spleen_56", "volume": "spleen_56.nii.gz", "text": "MSD Task09 spleen_56"}

A volume value can also name a folder holding one DICOM series. Give case files and folders pseudonymous names, because file and folder names reach the browser.

The config

The config below is complete and passes python -m potato.validate_cli config.yaml --strict on 2.12.1 and on 2.13.1:

annotation_task_name: Spleen segmentation, three readers
task_dir: .
output_annotation_dir: annotation_output/
media_directory: media
data_files:
  - data/cases.jsonl
item_properties:
  id_key: id
  text_key: text
user_config:
  allow_all_users: true
  users: []
 
# Every reader segments every case
num_annotators_per_item: 3
 
annotation_schemes:
  - annotation_type: volume_annotation
    name: spleen
    description: Paint the spleen on every axial slice that shows the spleen
    source_field: volume
    labels:
      - {name: Spleen, color: "#e6553a", key_value: "1"}
    # Boxes are saved as separate shapes and never become mask voxels
    tools: [brush, eraser]
    brush_size: 6
    # A NIfTI file has no modality header, so the CT windows are named here
    window_presets:
      - {name: Abdomen, center: 50, width: 400}
      - {name: Liver, center: 30, width: 150}
    default_window: Abdomen

Three settings need a reason:

• tools: [brush, eraser]. The viewer also has box and polygon tools. A box in the volume viewer is saved as a separate shape, never becomes mask voxels, and reaches neither the Dice score nor the NIfTI label volume. From Potato 2.13.1, closing a polygon paints the polygon's region into the mask on that slice, and Potato 2.12.1 saved a polygon as a shape. The readers in this run used the brush. The Navigate tool is always added.
• brush_size: 6. The brush radius is in voxels, with a default of 2 and a slider up to 15. At these scans' 0.7 to 0.8 mm pixels, a radius of 6 paints a disk about 9 mm across.
• window_presets. Potato offers its built-in CT windows only for a DICOM series whose header says CT. A NIfTI file has no such header, and without presets the viewer opens spleen_3 in an automatic window centered at -380 with a width of 1,287. The automatic window spreads 1,287 Hounsfield units over the gray scale and leaves little contrast between the spleen and the neighboring organs. The two presets are Radiopaedia's abdominal soft tissue and liver windows.

Reading the scans

Start the server in the project folder:

potato start config.yaml -p 8000

Each reader opens the server's address, registers with a username and a password, and gets all three cases, because num_annotators_per_item: 3 gives every case to three readers. The viewer shows axial, coronal, and sagittal views in radiological convention, with the patient's right on the viewer's left. The window menu opens on Abdomen (50/400) and also lists Liver (30/150) and the automatic window.

Press b for the brush and 1 for the Spleen label, then paint the spleen slice by slice. A stroke on the coronal or sagittal view writes into every axial slice the stroke crosses, and the eraser (e) removes voxels. Painting every slice is not required. Paint the first and last axial slice of the spleen and a few slices between the first and last where the outline changes, then press i to fill the slices between with Potato's shape-based interpolation. The filled slices are ordinary mask voxels, and Undo removes the filled slices. Clicking Next saves the case.

Agree on written rules before reading starts, such as whether vessels at the splenic hilum belong to the spleen, and give every reader the same copy. The agreement numbers below measure those rules as much as the readers' care. A rule added halfway through leaves the earlier cases read under a different rule.

The test run behind this page

The numbers on this page come from three scripted readers who worked through the running project and saved each case through the same requests the viewer sends. Each reader's mask came from the MSD expert mask, changed in one of three ways:

• ana painted every slice, one voxel inside the expert outline.
• ben painted the first slice, the last slice, and every third slice in between exactly, then pressed Interpolate. On spleen_3, ben painted 6 of the 14 slices and interpolation filled the other 8.
• cai painted every slice three voxels outside the expert outline and missed the top and bottom slice of the spleen.

The scripted readers make the result checkable against the expert masks, and a person's errors are less regular than the three scripted patterns.

Agreement between readers in Potato

Potato's agreement report scores volume_annotation masks by 3D Dice for each pair of readers on each case and each label. A pair where neither reader painted a label is left out for that label rather than counted as agreeing. Open /admin in a browser, type the admin key from admin_api_key.txt in the project folder, and then open /admin/iaa. For the test run, the report listed the three cases with all three readers each and gave:

mean_dice                0.887
dice_by_label.Spleen     0.887
n_pairs_by_label.Spleen  9

The 0.887 is the mean over 9 pairs, three per case. From Potato 2.13.1, the report also runs STAPLE on each case with three or more readers and averages each reader's rates over the cases. Run on this page's saved masks, the 2.13.1 report added:

staple.Spleen.sensitivity  ana 0.959  ben 0.939  cai 0.935
staple.Spleen.specificity  ana 1.000  ben 0.998  cai 0.944

The sensitivities equal the means of the SimpleITK figures further down. On a volume of more than two million voxels, Potato runs STAPLE inside the box around every reader's voxels, so the specificity counts only the voxels in that box and runs lower than the whole-scan figures from SimpleITK. The report gives no per-pair Dice. Computed from the exported masks, the pairs ranged from 0.954 (ana and ben on spleen_3) to 0.785 (ben and cai on spleen_56), and every pair with cai scored lower than ana with ben on the same case. Dice is not corrected for chance. Potato's chance-corrected geometry agreement covers boxes and polygons on 2D images and does not score volumes. The reliability export covers choice and rating questions only, and on this config python -m potato.export -c config.yaml -f reliability -o export/ stops with "No question type with a reliability coefficient".

Exporting NIfTI label volumes

The volume_nifti export writes one label volume per reader and case:

python -m potato.export -c config.yaml -f volume_nifti -o export/

export/labels__spleen.json
export/ana/spleen_3__spleen.nii.gz
export/ana/spleen_8__spleen.nii.gz
export/ana/spleen_56__spleen.nii.gz
export/ben/...
export/cai/...

Each .nii.gz file is a uint8 volume with 0 for background and 1 for Spleen, in the source scan's voxel order and affine. In the test run, every exported volume had the same shape and affine as the scan the volume came from. The masks therefore line up voxel for voxel with the scans. labels__spleen.json maps the values to label names. With more than one label in a scheme, a voxel that two labels share goes to the later label in the scheme's list.

For a DICOM series, -f dicom_seg writes a DICOM Segmentation that references the source series. On the release's CT phantom example, the export wrote one binary SEG with the annotator's username in the content description. The DICOM standard requires a SEG to copy the patient and study modules of the source series, and the phantom's SEG kept the phantom's patient name. Treat a SEG export as identifiable for that reason. On the NIfTI scans in this tutorial, dicom_seg writes no files and prints a warning to use volume_nifti.

STAPLE consensus and Dice for each reader

The script below runs SimpleITK's STAPLE filter on the three readers' masks for each case, keeps voxels with a consensus probability of 0.5 or more, and writes the consensus beside the readers' exports. The script then prints STAPLE's sensitivity and specificity for each reader, and each reader's Dice against the consensus. Install SimpleITK with pip install SimpleITK (the test run used 2.5.6), and run the script in the project folder:

from pathlib import Path
 
import SimpleITK as sitk
 
EXPORT = Path("export")
READERS = sorted(p.name for p in EXPORT.iterdir() if p.is_dir() and p.name != "staple")
CASES = sorted({f.name.split("__")[0] for r in READERS for f in (EXPORT / r).glob("*__spleen.nii.gz")})
(EXPORT / "staple").mkdir(exist_ok=True)
 
 
def dice(a, b):
    overlap = sitk.LabelOverlapMeasuresImageFilter()
    overlap.Execute(a, b)
    return overlap.GetDiceCoefficient()
 
 
print(f"{'case':10} {'reader':6} {'sens':>6} {'spec':>8} {'dice':>6}")
summary = {}
for case in CASES:
    masks = [sitk.ReadImage(str(EXPORT / r / f"{case}__spleen.nii.gz"), sitk.sitkUInt8) for r in READERS]
    staple = sitk.STAPLEImageFilter()
    staple.SetForegroundValue(1)
    probability = staple.Execute(masks)
    consensus = sitk.Cast(probability >= 0.5, sitk.sitkUInt8)
    consensus.CopyInformation(masks[0])
    sitk.WriteImage(consensus, str(EXPORT / "staple" / f"{case}__spleen.nii.gz"))
    for reader, mask, sens, spec in zip(READERS, masks, staple.GetSensitivity(), staple.GetSpecificity()):
        d = dice(mask, consensus)
        summary.setdefault(reader, []).append(d)
        print(f"{case:10} {reader:6} {sens:6.3f} {spec:8.5f} {d:6.3f}")
print("mean Dice against STAPLE:", {r: round(sum(v) / len(v), 3) for r, v in summary.items()})

case       reader   sens     spec   dice
spleen_3   ana     0.966  1.00000  0.984
spleen_3   ben     0.948  0.99997  0.971
spleen_3   cai     0.960  0.99907  0.901
spleen_56  ana     0.951  1.00000  0.980
spleen_56  ben     0.921  0.99998  0.954
spleen_56  cai     0.860  0.99947  0.834
spleen_8   ana     0.960  1.00000  0.981
spleen_8   ben     0.948  0.99993  0.971
spleen_8   cai     0.985  0.99795  0.921
mean Dice against STAPLE: {'ana': 0.982, 'ben': 0.965, 'cai': 0.885}

Read sensitivity and Dice first. Background voxels outnumber spleen voxels by 74 to 1 on spleen_8 and by 378 to 1 on spleen_56. Every reader's specificity therefore rounds to 0.998 or higher, and the small differences still rank the readers by how far each reader painted outside the spleen. Cai painted outside the outline and had the lowest specificity on all three scans. Cai's sensitivity on spleen_56 was 0.860, the lowest in the table, because the two slices cai missed were a larger share of a spleen that spans only 9 slices. Ben's interpolated masks scored between ana and cai on every case.

On these three scans, the STAPLE consensus was the same set of voxels as a two-of-three majority vote. With three readers and no reader far worse than the others, STAPLE's weights did not change any voxel's outcome. STAPLE adds the per-reader sensitivity and specificity in that case, and the consensus itself can differ from a vote with more readers or when careless readers outvote a careful reader.

Checking the consensus against the expert masks

The MSD expert masks make the test run checkable. Each reader's mask and the STAPLE consensus were scored by Dice against the expert mask:

The consensus scored higher against the expert mask than any single reader on every scan. Ana drew one voxel inside the outline and cai drew three voxels outside, and the two-of-three consensus fell between the two outlines. Scripted errors in opposite directions favor a consensus, and a study of real readers needs its own expert reference before claiming the same result.

Other tools for CT segmentation

3D Slicer is a free desktop application for detailed segmentation by one person. The Segment Editor's Fill between slices and Grow from seeds effects complete a segmentation from a few painted slices or seeds, and the SlicerRT extension's Segment Comparison module computes Dice and Hausdorff distances between two segmentations. ITK-SNAP is a free desktop application for 3D segmentation with semi-automatic active contours, and version 4.4.0 added AI features. MONAI Label is an Apache-2.0 server for AI-assisted labeling and active learning, used from viewers such as 3D Slicer, MITK, and OHIF. OHIF Viewer is an MIT-licensed web viewer that reads studies from a DICOMweb archive and supports segmentation workflows.

Potato fits when the same cases go to several readers through a browser on a server you run, and when the study needs the readers' agreement. Readers who prefer Slicer's or ITK-SNAP's segmentation tools can open Potato's NIfTI exports there. Medical Image and Clinical Text Annotation Tools Compared sets these tools beside Potato in more detail.

Limits of Potato for CT segmentation

• No consensus mask in Potato. Potato 2.13.1 reports pairwise Dice and STAPLE's rates for each reader, and the consensus mask comes from SimpleITK on the exported masks.
• No per-pair Dice. /admin/iaa gives one mean Dice per label over all pairs. Compute per-pair numbers from the exports, as above.
• No AI tool in the volume viewer. Masks come from the brush, the eraser, closed polygons, interpolation, and RTSTRUCT contours named on an item. MONAI Label, ITK-SNAP, and OHIF offer model-assisted 3D segmentation.
• No 4D volumes, fused PET/CT, or oblique reformats. A NIfTI file with more than one volume is refused, and an oblique acquisition is shown on the nearest orthogonal grid.
• No de-identification. Potato hides DICOM header fields from readers unless the config lists them, and leaves files, file names, and pixels as uploaded.

Where to go next

• Healthcare and clinical annotation, for Potato's clinical text and imaging features in one place
• Medical Image Annotation with Potato, for DICOM series, RTSTRUCT starting masks, and DICOM header filtering
• Medical Image and Clinical Text Annotation Tools Compared, for 3D Slicer, ITK-SNAP, MONAI Label, OHIF, and paid platforms beside Potato
• Export Formats, for volume_nifti, dicom_seg, and the other formats
• Task Assignment, for how Potato distributes cases among readers
• Air-Gapped Deployment, for running Potato on a network with no route out

References

1. Antonelli, M., Reinke, A., Bakas, S., Farahani, K., Kopp-Schneider, A., Landman, B. A., Litjens, G., Menze, B., Ronneberger, O., Summers, R. M., van Ginneken, B., Bilello, M., Bilic, P., Christ, P. F., Do, R. K. G., Gollub, M. J., Heckers, S. H., Huisman, H., Jarnagin, W. R., et al. (2022). The Medical Segmentation Decathlon. Nature Communications, 13, 4128. https://doi.org/10.1038/s41467-022-30695-9
2. Warfield, S. K., Zou, K. H., & Wells, W. M. (2004). Simultaneous truth and performance level estimation (STAPLE): An algorithm for the validation of image segmentation. IEEE Transactions on Medical Imaging, 23(7), 903-921. https://doi.org/10.1109/TMI.2004.828354