Calculation of GLSZM Failed.

[Haoxin22]

Hi,

I'm currently using the Pyradiomics v3 and when I try to calculate voxel-based features, there is a ""Calculation of GLSZM Failed."" error happened and seems like the error directs to the compiled C code which I cannot trace for.
###Image Format: .nrrd.
###Image protocol: multiple 2D MRI images, stack as a 3D image. E.g. 320x320x12 (horizontal axis pixel spacing is different from vertical axis pixel spacing and I'm not sure if that's the reason why it influences the result or not )
###Images and Image masks are valid. Image masks are also stacked 2D binary mask images, e.g. 320x320x12
###Log file:

DEBUG:radiomics.featureextractor: Parameters parsed, input is valid.
DEBUG:radiomics.featureextractor: Applying settings
DEBUG:radiomics.featureextractor: Enabled image types: {'Original': {}}
DEBUG:radiomics.featureextractor: Enabled features: {'glszm': None, 'gldm': None}
DEBUG:radiomics.featureextractor: Settings: {'binWidth': 25, 'label': 1, 'interpolator': 'sitkBSpline', 'resampledPixelSpacing': None, 'weightingNorm': None}
INFO:radiomics.featureextractor: Starting voxel based extraction
INFO:radiomics.featureextractor: Calculating features with label: 1
DEBUG:radiomics.featureextractor: Enabled images types: {'Original': {}}
DEBUG:radiomics.featureextractor: Enabled features: {'glszm': None, 'gldm': None}
DEBUG:radiomics.featureextractor: Current settings: {'minimumROIDimensions': 2, 'minimumROISize': None, 'normalize': False, 'normalizeScale': 1, 'removeOutliers': None, 'resampledPixelSpacing': None, 'interpolator': 'sitkBSpline', 'preCrop': False, 'padDistance': 5, 'distances': [1], 'force2D': False, 'force2Ddimension': 0, 'resegmentRange': None, 'label': 1, 'additionalInfo': True, 'binWidth': 25, 'weightingNorm': None, 'voxelBased': True}
INFO:radiomics.featureextractor: Loading image and mask
DEBUG:radiomics.imageoperations: Force casting mask to UInt32 to ensure correct datatype.
DEBUG:radiomics.imageoperations: Checking mask with label 1
DEBUG:radiomics.imageoperations: Calculating bounding box
DEBUG:radiomics.imageoperations: Checking minimum number of dimensions requirements (2)
DEBUG:radiomics.featureextractor: Image and Mask loaded and valid, starting extraction
DEBUG:radiomics.featureextractor: Creating image type iterator
INFO:radiomics.featureextractor: Adding image type "Original" with custom settings: {}
DEBUG:radiomics.featureextractor: Extracting features
DEBUG:radiomics.imageoperations: Yielding original image
INFO:radiomics.featureextractor: Calculating features for original image
DEBUG:radiomics.imageoperations: Cropping to size [72 76 12]
INFO:radiomics.featureextractor: Computing glszm
DEBUG:radiomics.glszm: Initializing feature class
DEBUG:radiomics.imageoperations: Discretizing gray levels inside ROI
DEBUG:radiomics.imageoperations: Calculated 34 bins for bin width 25 with edges: [ 25 50 75 100 125 150 175 200 225 250 275 300 325 350 375 400 425 450
475 500 525 550 575 600 625 650 675 700 725 750 775 800 825 850 875])
DEBUG:radiomics.glszm: Calculating voxel batch no. 1/1
DEBUG:radiomics.glszm: Calculating GLSZM matrix in C

###Error Message:
Traceback (most recent call last):

File "mypath/making_radiomics_dataset_with_biomarker.py", line 188, in
main()

File "mypath/making_radiomics_dataset_with_biomarker.py", line 133, in main
result_t2 = extractor.execute(curr_t2_path, curr_roi_path, voxelBased=True)

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/featureextractor.py", line 332, in execute
featureVector.update(self.computeFeatures(inputImage, inputMask, imageTypeName, **inputKwargs))

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/featureextractor.py", line 530, in computeFeatures
for (featureName, featureValue) in six.iteritems(featureClass.execute()):

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/base.py", line 183, in execute
self._calculateVoxels()

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/base.py", line 208, in _calculateVoxels
for success, featureName, featureValue in self._calculateFeatures(voxelCoords):

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/base.py", line 229, in _calculateFeatures
self._initCalculation(voxelCoordinates)

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/glszm.py", line 65, in _initCalculation
self.P_glszm = self._calculateMatrix(voxelCoordinates)

File "/home/haoxin/anaconda3/lib/python3.6/site-packages/pyradiomics-0+unknown-py3.6-linux-x86_64.egg/radiomics/glszm.py", line 93, in _calculateMatrix
P_glszm = cMatrices.calculate_glszm(*matrix_args) # shape (Nvox, Ng, Ns)

IndexError: Calculation of GLSZM Failed

Would you mind helping me? The shape, firstorder, GLCM, GLRLM are calculated just OK. I cannot figure out why the calculation of GLSZM failed since seems like the error happens inside cMatrices and it directs to _cmatrices which is a private function that I cannot get access to.

Thanks

[JoostJM]

That is indeed an error from inside the C extension. Can you share the image, mask and parameters that caused this error? That way, I may be able to track down this bug. Does it happen for just 1 case or multiple cases in your dataset?

[Haoxin22]

Thank you so much for your reply. Sure, of course. Please see as following:

1. I've tried different patients, and the error happened for every patient in my dataset.

2. Image&mask settings: The original Image is in dicom format and all images are 320x320x12 shape. I used pydicom read the image into numpy array and use nrrd.write to write the 3-dim nparray into .nrrd as 3D image output. The mask with the same dimension and the same process, the only difference is that the mask contains only 0(not RoI)/1(RoI).
   My possible question regarding this part: (1) is that possible the transformation from dicom->nrrd lose some info? (2) Is that possible that the datatype matters? I didn't find the datatype requirements though.

3. Image examples & mask examples & params settings:

https://drive.google.com/drive/folders/1XFTM5KoeKEqMPlVF83sFyqpo_v2IzZj_?usp=sharing

The link contains 2 image&mask examples, and 2 settings of parameters (both do not work). And also the image datatype print (the same as mask img datatype). I didn't do any preprocessing to the image, just the original image from dicom file and stored in nrrd format. If I set the "voxelBased = False", everything goes well. But if I set the "voxelBased = True", the mentioned error happened.

Hope that helps! If you need the actual data example but not screen shot, please let me know.

Two more further question: for the 3D MR data like I have, 320x320x12, whose resolution in x-y direction is higher than z direction,

Q1: if I'd like to get the radiomics features for the whole 3D MR image, then should I set the voxelBased = True or False? For example, I only want to get one firstorder_10Percentile value for the whole 3D data but not 12 different firstorder_10Percentile values. Not sure what's the difference for voxelBased = True or False. I thought for 3D it should be True since we are facing 3D image and the unit is voxel. If set to False, will the radiomics features be calculated for each slice and then be averaged?

Q2: Why the returned value when voxelBased = True is a feature map? I'm not quiet understand how to compute radiomics features for EACH VOXEL. For example, how to calculate firstorder_10Percentile for the voxel [200,200,6]?

Thank you so much!

[jvanlunenburg]

It seems you haven't fully understood the voxelBased setting. By definition, voxelBased is a feature map (aka parametric map). If false, it returns a "segment-based" scalar per feature, according to the documentation.

if I'd like to get the radiomics features for the whole 3D MR image, then should I set the voxelBased = True or False?

So you should set it to false.

I'm not quiet understand how to compute radiomics features for EACH VOXEL.

Then you should set it to true. What it does is use the local pixel neighbourhood (e.g. 3x3x3) to get per-pixel textures.

[Haoxin22]

Thank you so much. One more question, how does it calculate the radiomics for a 3D image if I set voxelBased = False? Does it calculate 2D features for each slice and then average them?

[jvanlunenburg]

Thank you so much. One more question, how does it calculate the radiomics for a 3D image if I set voxelBased = False? Does it calculate 2D features for each slice and then average them?

That depends on the feature category. The IBSI defines several possible types of "aggregation logic", the type you mentioned is one of them. If we google "pyradiomics GLCM" we end up on this page. The relevant part is:

By default, the value of a feature is calculated on the GLCM for each angle separately, after which the mean of these values is returned.

So to paraphrase based on my understanding, by default it computes a set of 3D GLCMs corresponding to a set of directions, then gets a feature per GLCM, then averages those to return a rotationally invariant scalar feature.

[JoostJM]

@jvanlunenburg correct. Though defined in IBSI, PyRadiomics does not support calculating features per 2D slice and then averaging them, as this would give unreasonable emphasis on voxels located on slices with only few voxels segmented. PyRadiomics does support "force2D", even with 3D input, but in that case, texture matrices are still calculated based on the 3D mask, across all slices, the difference is that directions that define an offset moving between slices are ignored.

E.g. GLCM 3D mask: 13 matrices for 13 directions, each defining a GLCM based on all segmented voxels
GLCM 3D mask, force2D=True: 4 matrices for 4 directions moving "in-plane", each defining a GLCM based on all segmented voxels.
Features are calculated on the matrices and then averaged across the directions. Optionally, matrices can be merged prior to feature calculation, this is done by specifying a valid value for weightingNorm.

[JoostJM]

@Haoxin22, can you please share the NRRD files of the image and mask for me to track down your issue?

[fer-ferreira]

Hi! I'm having the same issue...

T02C04_flair.nii.gz
T02C04_mask.nii.gz

[JoostJM]

Fixed by #635