CLIJ2 overlapMatrix and JaccardIndexMatrix XY coordinate issues

Hi,

I found an issue with two of CLIJ2’s overlap commands, namely CLIJ2_generateBinaryOverlapMatrix and CLIJ2_generateJaccardIndexMatrix both yield incorrect results for me WRT indexes and XY coordinates. I’ve written up a little script and included some test data to illustrate the problem. Just change the directory and it should run.

CLIJ_Overlap_Test.ijm (5.4 KB)
CLIJ_Overlap_Test_Data_8Bit.tif (616.7 KB)

Here’s a brief description of the issue, which I’ve also written out in the macro:
Assuming the setup is as follows:

Ext.CLIJ2_generateBinaryOverlapMatrix(class1_labelmap, class2_labelmap, class1_class2_overlap_matrix);
Ext.CLIJ2_generateJaccardIndexMatrix(class1_labelmap,  class2_labelmap,  class1_class2_jaccard_index_matrix);

generateBinaryOverlapMatrix adds not only the correct results to the output matrix, but also a mirrored&flipped version of these results. In other words, if class1 object4 overlaps with class2 object 2, the matrix has a value of 1 at (4,2), which is the correct result, but also at (2,4), an incorrect, mirrored result. This issue can be observed by adding a scalar to the one of the lablemaps (here class2), which shifts the indexes to more positive values of the corresponding axis. You can get the correct matrix by adding the number_of_class1 to the class2_labelmap, then cropping the result.

generateJaccardIndexMatrix reverses the X and Y coordinates for the object indexes, i.e. it adds the X results to Y coordinate and vice versa. For example, if class1 object4 overlaps with class2 object 2, the matrix should have a value non-zero value at (4,2), but instead has it at (2,4). Again, this becomes apparent when adding a scalar to class2_labelmap to shift the indexes. The result should be shifted in the Y direction, with the elongation of this axis. Instead, the Y axis is elongated as expected, but the results are shifted by a corresponding amount along the X axis.

I’ve found a way to take care of the generateBinaryOverlapMatrix issue, but not the generateJaccardIndexMatrix one. Either way, wanted to make the devs aware of it in case others run into the same problems. Also interested if others have had the same issues and/or know some workarounds/solutions. Any help is appreciated!

Patrick

Some other (possibly) relevant info:
I’m running the most up-to-date version of CLIJ/CLIJ2 on imageJ v1.53c.
Windows 7
GPU: Intel® HD Graphics 530
Memory in GB: 1.2672
OpenCL version: 2

// Object overlap on GPU
// 	generateBinaryOverlapMatrix and generateJaccardIndexMatrix both yield problematic results. 
// 		> 	generateBinaryOverlapMatrix adds not only the correct results to the output matrix, but also a mirrored&flipped version of these results.
//			In other words, if class1 object4 overlaps with class2 object 2, the matrix has a value of 1 at (4,2), which is the correct result, but  
//			also at (2,4), an incorrect, mirrored result. This issue can be observed by adding a scalar to the one of the lablemaps (here class2),
//			which shifts the indexes to more positive values of the corresponding axis. The correct result is achieved by adding the number_of_class1
//			to the class2_labelmap, then cropping the result. 
//		>	generateJaccardIndexMatrix reverses the X and Y coordinates for the object indexes, i.e. it adds the X results to Y coordinate and vice versa. 
//			Again, this becomes apparent when adding a scalar to class2 to shift the indexes. The result should be shifted in the Y direction, with the elongation of this
//			axis. Instead, the Y axis is elongated as expected, but the results are shifted by a corresponding amount along the X axis.



	requires("1.53c");
	run("CLIJ2 Macro Extensions", "cl_device=");
	Ext.CLIJ2_getGPUProperties(gpu, memory, opencl_version);
	print("GPU: " + gpu);
	print("Memory in GB: " + (memory / 1024 / 1024 / 1024) );
	print("OpenCL version: " + opencl_version);	
	
// Clean up at the beginning
	close("*");
	Ext.CLIJ2_clear();
	run("Clear Results");
	saveSettings();

	print("");
	
// Loading images, extracting dimensional data, and pushing them into GPU memory
	dir = // "directory/of/test/image/";
	image = "CLIJ_Overlap_Test_Data_8Bit.tif";
	name = substring(image, 0, indexOf(image, "."));
	path = dir+image;
	run("Bio-Formats Importer", "open=["+path+"] autoscale color_mode=Default concatenate_series open_all_series rois_import=[ROI manager] view=Hyperstack stack_order=XYCZT");
	rename("CLIJ_Overlap_Test_Data_8Bit");
	getDimensions(width, height, channels, slices, frames);
	imageWidth = width;
	imageHeight = height;
	if (nSlices == 1){
		nSpatialDimensions = 2;
	}else{
		nSpatialDimensions =  3;
	}	
	run("Split Channels");
	selectWindow("C1-"+name);
	class1_mask = getTitle();
	Ext.CLIJ2_push(class1_mask);		
	selectWindow("C2-"+name);
	class2_mask = getTitle();
	Ext.CLIJ2_push(class2_mask);

	close("*");

// Label connected points (generate superpixels), find centroids to get object counts
	
	class1_centroids = "class1_centroids";
	class1_labelmap = "class1_labelmap";
	Ext.CLIJ2_connectedComponentsLabelingBox(class1_mask, class1_labelmap);
	Ext.CLIJ2_centroidsOfLabels(class1_labelmap, class1_centroids);
	Ext.CLIJ2_getSize(class1_centroids);
	/// object count	
	number_of_class1 = getResult("Width", nResults() - 1);
	nFrames = getResult("Depth", nResults() - 1);
	IJ.log("number of dimensions: " + nSpatialDimensions);
	IJ.log("number of frames: " + nFrames);
	print("");
	IJ.log("number of class1: " + number_of_class1);
	/// pull labelmap
	Ext.CLIJ2_pull(class1_labelmap);
	run("glasbey on dark");
	setMinAndMax(0, number_of_class1);
	if(nSpatialDimensions == 3){Stack.setSlice(nSlices/3);}

	// class2 
	/// labelmap & centroid
	class2_centroids = "class2_centroids";
	class2_labelmap = "class2_labelmap";
	Ext.CLIJ2_connectedComponentsLabelingBox(class2_mask, class2_labelmap);
	/// object count	
	Ext.CLIJ2_centroidsOfLabels(class2_labelmap, class2_centroids);
	Ext.CLIJ2_getSize(class2_centroids);
	number_of_class2 = getResult("Width", nResults() - 1);
	IJ.log("number of class2 : " + number_of_class2);
	/// pull labelmap
	Ext.CLIJ2_pull(class2_labelmap);
	run("glasbey on dark");
	setMinAndMax(0, number_of_class2);
	if(nSpatialDimensions == 3){Stack.setSlice(nSlices/3);}

	Ext.CLIJ2_release(class1_centroids);
	Ext.CLIJ2_release(class2_centroids);
	
	run("Clear Results");

// create an overlap matrix

	class1_class2_overlap_matrix = "class1_class2_overlap_matrix";
	Ext.CLIJ2_generateBinaryOverlapMatrix(class1_labelmap, class2_labelmap, class1_class2_overlap_matrix);
	Ext.CLIJ2_pull(class1_class2_overlap_matrix);
	run("HiLo");
	class1_class2_jaccard_index_matrix = "class1_class2_jaccard_index_matrix";
	Ext.CLIJ2_generateJaccardIndexMatrix(class1_labelmap,  class2_labelmap,  class1_class2_jaccard_index_matrix);
	Ext.CLIJ2_pull(class1_class2_jaccard_index_matrix);
	run("HiLo");

	scalarArray = newArray(Math.floor(number_of_class1/2), number_of_class1);
	for(i=0; i<scalarArray.length; i++){
		if(i>0){Ext.CLIJ2_release(class2_labelmap_adj);}
		scalar = scalarArray[i];
		Ext.CLIJ2_addImageAndScalar(class2_labelmap, class2_labelmap_adj, scalar);
		class1_class2_overlap_matrix_adj = "class1_class2plus"+scalar+"_overlap_matrix";
		Ext.CLIJ2_generateBinaryOverlapMatrix(class1_labelmap, class2_labelmap_adj, class1_class2_overlap_matrix_adj);
		Ext.CLIJ2_pull(class1_class2_overlap_matrix_adj);
		run("HiLo");
		if(i==scalarArray.length-1){
			expected_result = "expected_result";
			Ext.CLIJ2_crop2D(class1_class2_overlap_matrix_adj, expected_result, 0, number_of_class2-1, number_of_class1+1, number_of_class2+1);
			Ext.CLIJ2_pull(expected_result);
			run("HiLo");
		}
		class1_class2_jaccard_index_matrix_adj = "class1_class2plus"+scalar+"_jaccard_index_matrix"+scalar;
		Ext.CLIJ2_generateJaccardIndexMatrix(class1_labelmap,  class2_labelmap_adj,  class1_class2_jaccard_index_matrix_adj);
		Ext.CLIJ2_pull(class1_class2_jaccard_index_matrix_adj);
		run("HiLo");
	}


	Ext.CLIJ2_clear();
1 Like

Hey Patrick @pdd2110,

thanks for reporting. This is obviously a bug!

First of all the binary overlap matrix: This little code snippet

run("CLIJ2 Macro Extensions", "cl_device=");
Ext.CLIJ2_clear();
print("\\Clear");

labels_1 = "1 0 2 2 3 0 0 4";
labels_2 = "1 0 5 2 2 3 0 4";

// --------------------------------------
// push data to GPU
width = lengthOf(labels_1);
height = 1;
depth = 1;

Ext.CLIJ2_pushString(label_map_1, labels_1, width, height, depth);
Ext.CLIJ2_pushString(label_map_2, labels_2, width, height, depth);

// generate binary overlap matrix
Ext.CLIJ2_generateBinaryOverlapMatrix(label_map_1, label_map_2, binary_overlap_matrix);

// show the matrix in the log window
Ext.CLIJ2_print(binary_overlap_matrix);

outputs this binary overlap matrix:

1.0 0.0 0.0 1.0 0.0
0.0 1.0 0.0 0.0 0.0
0.0 0.0 1.0 1.0 0.0
1.0 0.0 1.0 0.0 0.0
0.0 0.0 0.0 0.0 1.0
0.0 0.0 1.0 0.0 0.0

As you pointed out correctly, it contains a mirrored copy of the actually correct matrix. I just fixed this bug and it now outputs:

1.0 0.0 0.0 0.0 0.0
0.0 1.0 0.0 0.0 0.0
0.0 0.0 1.0 1.0 0.0
1.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 1.0
0.0 0.0 1.0 0.0 0.0

Regarding the Jaccard index matrix. A similar example script:

run("CLIJ2 Macro Extensions", "cl_device=");
Ext.CLIJ2_clear();
print("\\Clear");

labels_1 = "1 0 2 2 3 0 0 4";
labels_2 = "1 0 5 2 2 3 0 4";

// --------------------------------------
// push data to GPU
width = lengthOf(labels_1);
height = 1;
depth = 1;

Ext.CLIJ2_pushString(label_map_1, labels_1, width, height, depth);
Ext.CLIJ2_pushString(label_map_2, labels_2, width, height, depth);

// generate binary overlap matrix
Ext.CLIJ2_generateJaccardIndexMatrix(label_map_1, label_map_2, overlap_matrix);

// show the matrix in the results table window
run("Clear Results");
Ext.CLIJ2_pullToResultsTable(overlap_matrix);

outputs this table:

image

It takes some time to see it, but this matrix is x-y-transposed to what it is supposed to be. I just fixed that bug as well and the matrix looks now like this:

image

Let’s quickly go through it to make sure I got it right this time:

  • Label 0/0 (Column X0, Row 1): There are three pixels in question (0 in label_1 or 0 in label_2), two of them match, J=0.66.
  • Label 1/1 (Column X1, Row 2): There is one pixel in question and it matches, J=1
  • Label 2/2 (Column X2, Row 3): There are three pixels in question (2 in label_1 or 2 in label_2), one matches, J=0.33
  • Label 2/3 (Column X2, Row 4): There are three pixels in question, none matches, J=0
  • Label 3/2 (Column X3, Row 3): There are 2 pixels in question, one matches, J=0.5
  • Label 3/3 (Column X3, Row 4): There are 2 pixels in question, none matches, J=0
  • Label 4/4 (Column X4, Row 5); There is one pixel in question, it matches, J=0
  • Label 2/5 (Column X2, Row 6): There are two pixels in question, one matches, J=0.5

Let me know if you see any other discrepancy.

Thanks a lot for reporting these bugs! You made CLIJ better today. I’ll send around an update asap and you should be able to try it out by tomorrow. :slight_smile:

Thanks again!

Cheers,
Robert

1 Like

Thanks a ton Robert, glad I could help! I’ll let you know if I run into anything else. Really enjoying the ease (and power) of CLIJ, and looking forward to the update! Thanks for all of the hard work!

1 Like

Ok, the update is out. Let me know if it works now for your use-case!

Thanks again!

Cheers,
Robert

Seeing your response here a bit late. I updated last night and it’s been working perfectly!

Thanks for being so responsive!
Patrick

1 Like