Colour Gradient Based on Pixel Density of Overlaid Segmentations?

Support
,
1

Hello,

I am new to 3D slicer, but I have a series of MRIs which I have overlaid and segmented. I would like to now take those segmentations and apply a filter to visualize their overlap. I have seen this done with two segmentations, but as I have ~40 I would like to use a colour scale to express the relative amount of overlap.

Any suggestions would be appreciated!

Thanks, Greg

2

A simple way of dealing with this is to write a short Python script, which does the following:

  • create an “accumulation volume” (segments from all images will be added to it) fill with all zeros, get it as a numpy array (using slicer.util.arrayFromVolume())
  • load each segmentation, export the segment to a numpy array (using slicer.util.arrayFromSegmentBinaryLabelmap), you threshold this numpy array (set it to 1 where it has non-zero values), and add the values to the accumulation volume’s numpy array
  • update the accumulation volume node from its numpy array (using slicer.util.arrayFromVolumeModified())
3

Hi Andras,

Thank you for your response. Sorry I am still a little confused, I believe I have followed your instructions (although I am very much still learning so the mistake is probably mine), but now I just see a complete overlay of all these segmentations. Is there a way to apply a colour gradient based on the relative overlap?

Thanks, Greg

4

Could you copy the otthon script that you created? We could tell you what to add or change to see the result as a scalar volume instead of a binary.

5

Hi Andras,

I am not terribly familiar with python (so after a few failed attempts) I used the method you mention under: Segment / BinaryLabelMap to numpy array.

Would this work for the problem?

6

You need to add the numpy arrays so that if a voxel is set to 1 in multiple segmentations then that voxel value will be larger than 1 (the value should be the number of segmentations where that voxel is non-zero).

7

Hi Andras,

I think I must be making some mistake, when I try to make an accumulation volume using the code above I get thrown:

File “”, line 1, in
File “/Applications/Slicer.app/Contents/bin/Python/slicer/util.py”, line 1476, in arrayFromVolume
vimage = volumeNode.GetImageData()
AttributeError: ‘str’ object has no attribute ‘GetImageData’

Greg

8

arrayFromVolumeNode takes a volume node object as input (not the volume’s name). You can get the volume node object from the volume’s name by calling:

volumeNode = slicer.util.getNode('MyVolumeName')
9

Hi Andras,

When entering this code I run into the following problem:
Traceback (most recent call last):
File “”, line 1, in
File “/Applications/Slicer.app/Contents/bin/Python/slicer/util.py”, line 1324, in getNode
raise MRMLNodeNotFoundException(“could not find nodes in the scene by name or id ‘%s’” % (pattern if (isinstance(pattern, str)) else “”))
slicer.util.MRMLNodeNotFoundException: could not find nodes in the scene by name or id ‘the_scene_I_inputted.mrml’.

Sorry for all of these small errors

10

Replace MyVolumeName with the actual name of your volume node. It seems that the string that the method got was the_scene_I_inputted.mrml which does not seem like a volume node name.

11

Hi Andras, sorry I should have clarified.

I replaced ‘My Volume Name’ with the volume of the scene I was working with. The ‘the_scene_I_inputted.mrml’ is meant to be the placeholder for the scene that I inputted under ‘My Volume’. Sorry for the confusion.

12

This code should work. It produces a volume named ExportVolume, which you can use in Volume Rendering to assign colors and opacities to show the density of original segments in each location.

# Export one volume that shows the sum of segments in a segmentation of the same shape

# Update these names using your scene

segmentationNodeName = '2'
masterVolumeName = '10 ciss3d_tra_iso'
outputVolumeName = 'ExportVolume'

# Prepare nodes that will be needed

segmentationNode = slicer.mrmlScene.GetFirstNodeByName(segmentationNodeName)
segLogic = slicer.modules.segmentations.logic()
segmentation = segmentationNode.GetSegmentation()
labelmapVolume = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLLabelMapVolumeNode")

# Make sure the output volume exists, and it has the same shape as the master volume for segmentation

outputVolume = slicer.mrmlScene.GetFirstNodeByName(outputVolumeName)
masterVolume = slicer.mrmlScene.GetFirstNodeByName(masterVolumeName)
if outputVolume is None:
  volumesLogic = slicer.modules.volumes.logic()
  outputVolume = volumesLogic.CloneVolume(slicer.mrmlScene, masterVolume, outputVolumeName)
  outputVolume.CreateDefaultDisplayNodes()

# Prepare output volume array and set voxel values to zero

exportArray = slicer.util.arrayFromVolume(outputVolume)
exportArray[:] = 0

# Add each segment to the output array one-by-one

for i in range(segmentation.GetNumberOfSegments()):
  segmentIds = vtk.vtkStringArray()
  segmentId = segmentation.GetNthSegmentID(i)
  segmentIds.InsertNextValue(segmentId)  
  segLogic.ExportSegmentsToLabelmapNode(segmentationNode, segmentIds, labelmapVolume, masterVolume)
  labelArray = slicer.util.arrayFromVolume(labelmapVolume)
  exportArray = exportArray + labelArray

# Update the output volume

slicer.util.updateVolumeFromArray(outputVolume, exportArray)