Convert 2D image to 1D Array in CLIJ2?

Hi!

I am trying to get statistics on the 10% brighter pixels of each slice of a stack. My idea was to create a 1D array for all the intensities in the image, and then find the 90th percentile using a custom function.

In the CLIJ2 plugin, I have found this super useful function “Convert Table Column to 2D image”, but not the other way around. In order to obtain a 1D array of intensities, I am using the CLIJ2 function “Convert 2D image to a table on GPU”, and then converting the table into a 1D array (see attached macro).

Unfortunately, this macro is pretty slow, so it takes ~8.5s to process the stack in attachment (50 slices), which is too slow, since I need to process hundreds of stacks of ~2000 slices. Do you have suggestions or ideas I can implement to speed up the process?

Thanks a lot!!
Emanuele

convert_slice_to_1d_ARRAY.ijm (1.3 KB)

myStack.tif (585.3 KB)

Hi!

I have actually found a workaround using getRawStatistics method, and this is much faster (it runs in ~0.9s, so x10 faster than the previous macro). I am attaching the macro.

However, it takes ~34s to run for the same stack with all 2000 slices. Is there an alternative way to further accelerate the process (maybe using some GPU acceleration)?

Thanks,
Emanuele

convert_slice_to_1d_ARRAY_hist.ijm (1.1 KB)

1 Like

Hey @eleo,

cool question. Thanks for asking! I just checked. If you change your workflow a bit, you can generate the histograms on the GPU and pull them back as a table before you do your detailed analysis. Not sure how long the second part takes, but collecting the histograms of 2000 image slices can be done in 2.5 seconds.

But one counter question: Wouldn’t it make sense to generate the histogram of the whole stack, find the intensity threshold corresponding to the 90th percentile and then analyse in all slices the distribution of the pixel intensities above that threshold? Not sure, this might depend on the physics of your experimental setup.

If you want to give it a try, this is a CLIJ2 macro:

open("C:/Users/rober/Downloads/myStack.tif");
start_time = getTime();

// init GPU
run("CLIJ2 Macro Extensions", "cl_device=[RTX]");
Ext.CLIJ2_clear();
Ext.CLIJ2_getGPUProperties(GPU_name, global_memory_in_bytes, OpenCL_version);
print("Running on " + GPU_name);

// push data to GPU
stack = getTitle();
Ext.CLIJ2_push(stack);

// do statistics on the stack
Ext.CLIJ2_getMaximumOfAllPixels(stack, maximum_of_all_pixels);

// allocate memory for the histogram and the histogram collection
number_of_bins = 256;
Ext.CLIJ2_create2D(histogram, number_of_bins, 1, 32);
Ext.CLIJ2_getDimensions(stack, width, height, depth);
Ext.CLIJ2_create2D(histogram_collection, number_of_bins, depth, 32);

// go through slices
for (z = 0; z < 2000; z++) {
	// get a slice from the stack
	Ext.CLIJ2_copySlice(stack, slice, z);

	// get the histogram from the slice
	Ext.CLIJ2_histogram(slice, histogram, 256, 0, maximum_of_all_pixels, false);

	// copy the histogram
	Ext.CLIJ2_paste2D(histogram, histogram_collection, 0, z);
}

// get all histograms in a table
run("Clear Results");
Ext.CLIJ2_pullToResultsTable(histogram_collection);

print("The whole process took " + (getTime() - start_time) + " ms");

If you execute it on Intel and NVidia GPUs, they almost have the same processing time. This is a hint that processing basically takes no time and push/pull are dominant:
image

Also, in order to exploit GPU-acceleration optimally, running whole workflow on the GPU makes a lot of sense. Single operations such as histogram generation might not pay off. Read more.

Last but not least: Thanks to your question, I found a bug in the histogram method of clij2. That’s why: Thanks for asking! And please update your Fiji, if you have clij2 installed already, before you run the script above.

Let me know if this helps!

Cheers,
Robert

1 Like

I just realised that you already found the right methods :slight_smile: The method I’m using “pullToResultsTable()” is just a new name for “convert2DImageToResultsTable()” - they do the same :wink:

A “pullArray” method is in the making, but not functional yet.

Thanks for testing CLIJ2!

1 Like

Hi Robert,

Thanks a lot for your reply and for the macro!!!

Your counter question is super interesting! I think I didn’t explain what is my problem in this case. I am trying to extract the timing information for when the signal from the cells get to its maximum (from the zero signal in the beginning of the stack to its maximum) without segmenting the cells. My idea was to consider the average of the top 10% of pixels of each slice as a measure of the intensity. This should be more robust than considering the maximum, and I should filter out variations in pixel intensities, but also this should be a more robust metric than considering the average intensity value of a rectangle in the center of the stack. This is why I was considering the percentile for each slice, instead of the percentile of the entire stack, since I thought that this would give me an abrupt increase, instead of a gradual increase from background to the maximum.

I implemented your code in my machine and it takes 3.8s for the 50 slices stack (myStack.tif). I learned quite a few new CLIJ2 methods! I have few questions on the CLIJ2 methods you used:

  1. What is the output/result CLIJ2_histogram? I think to understand you use a 2D image (1 x 256) to save the result of the histogram, and store it as the z-th vertical column of a 2D image (nSlices x 256) using CLIJ2_paste2D. Does it needs an image that is 1 x nBins to store properly the data? Is the column number the intensity and the value the ‘counts’ for that intensity?
  2. I would like to obtain a final table where each column represent the histogram count for each slice. I achieved it by rotating counterclockwise the image before applying the CLIJ2_pullToResultsTable method. Is it the right thing to do?
  3. I modified the number_of_bins = maximum_of_all_pixels+1; in order to have a count for all values in the image from 0 to the maximum in the stack. However, when I look at the result of the algorithm I get 0 count for all values beyond 124, while the maximum in the stack is 204. Did I do something wrong?

Here is a slightly modified version of the macro you sent me

open("C:/Users/rober/Downloads/myStack.tif");
start_time = getTime();

// init GPU
run("CLIJ2 Macro Extensions", "cl_device=[RTX]");
Ext.CLIJ2_clear();
Ext.CLIJ2_getGPUProperties(GPU_name, global_memory_in_bytes, OpenCL_version);
print("Running on " + GPU_name);

// push data to GPU
Ext.CLIJ2_push(stack);

// do statistics on the stack
Ext.CLIJ2_getMaximumOfAllPixels(stack, maximum_of_all_pixels);

// allocate memory for the histogram and the histogram collection
number_of_bins = maximum_of_all_pixels+1;
Ext.CLIJ2_create2D(histogram, number_of_bins, 1, 32);
Ext.CLIJ2_getDimensions(stack, width, height, depth);
Ext.CLIJ2_create2D(histogram_collection, number_of_bins, depth, 32);

// go through slices
for (z = 0; z < slices_num; z++) {
   // get a slice from the stack
   Ext.CLIJ2_copySlice(stack, slice, z);

   // get the histogram from the slice
   Ext.CLIJ2_histogram(slice, histogram, number_of_bins, 0, maximum_of_all_pixels, false);

   // copy the histogram
   Ext.CLIJ2_paste2D(histogram, histogram_collection, 0, z);
}
// rotate image in order to get the histogram of each slice as a column in the table
histogram_collection_rot = histogram_collection + "_ROTATE";
Ext.CLIJ2_rotateCounterClockwise(histogram_collection, histogram_collection_rot);

// get all histograms in a table
run("Clear Results");
Ext.CLIJ2_pullToResultsTable(histogram_collection_rot);

print("The whole process took " + (getTime() - start_time) + " ms");

Thanks a bunch!
Emanuele

PS: I’m glad that the question allowed to solve a bug! Helping without realizing it :slight_smile:

Good questions, especially as the documentation for the histogram method was missing. I just added it on the website. In short: A histogram has n bins and the corresponding image has width=n, height=1, depth=1. The intensity in that image corresponds to the number of pixels with gray value in a given range.

Sure! Maybe it’s better to use transposeXY instead of rotating because then the first histogram entry is on top of the table.

Two things: 1) This was apparently another bug in the histogram() function. Update your Fiji to get the recent version. 2) If you use rotateCounterClockwise, you rotate the end of the histogram to the top of the table. RotateClockwise or transposeXY solve this problem.

In order to test it, here comes the updated macro:

open("C:/Users/rober/Downloads/myStack.tif");
stack = getTitle();
start_time = getTime();

// init GPU
run("CLIJ2 Macro Extensions", "cl_device=[RTX]");
Ext.CLIJ2_clear();
Ext.CLIJ2_getGPUProperties(GPU_name, global_memory_in_bytes, OpenCL_version);
print("Running on " + GPU_name);

// push data to GPU
Ext.CLIJ2_push(stack);

// do statistics on the stack
Ext.CLIJ2_getMaximumOfAllPixels(stack, maximum_of_all_pixels);
print("max: " + maximum_of_all_pixels);

// allocate memory for the histogram and the histogram collection
number_of_bins = maximum_of_all_pixels+1;
Ext.CLIJ2_create2D(histogram, number_of_bins, 1, 32);
Ext.CLIJ2_getDimensions(stack, width, height, depth);
Ext.CLIJ2_create2D(histogram_collection, number_of_bins, depth, 32);

// go through slices
for (z = 0; z < depth; z++) {
    // get a slice from the stack
    Ext.CLIJ2_copySlice(stack, slice, z);

   // get the histogram from the slice
   Ext.CLIJ2_histogram(slice, histogram, number_of_bins, 0, maximum_of_all_pixels, false);

   // copy the histogram
   Ext.CLIJ2_paste2D(histogram, histogram_collection, 0, z);
}
// rotate image in order to get the histogram of each slice as a column in the table
Ext.CLIJ2_transposeXY(histogram_collection, histogram_collection_rot);

// get all histograms in a table
run("Clear Results");
Ext.CLIJ2_pullToResultsTable(histogram_collection_rot);

print("The whole process took " + (getTime() - start_time) + " ms");

Thanks again for testing! It’s a pleasure working with you!

Cheers,
Robert

Hi Robert,

Thanks again for all your explanations, for updating the library and for the updated macro!! I will test it and let you know about it!

Thanks for catching that mistake of rotating counter-clockwise. I realize the mistake I made there. :slight_smile: Concerning the histogram, if I want to count the number of occurrences of each intensity, I have simply to give the maximum intensity of the image as the number of bins, right?

Cheers,
Emanuele

1 Like

Hey @eleo,

You had it correct earlier:

But to clarift and remove any confusion, let’s try it out. I determine the histogram of a mini-image that has intensities between 0 and 5:

run("CLIJ2 Macro Extensions", "cl_device");


// create an image with intensities 0 .. 5
array = newArray(0,1,2,3,4,5,4,3,2,1);
Ext.CLIJ2_pushArray(image, array, 10, 1, 1);

// make a histogram with 5 bins
number_of_bins = 5;
Ext.CLIJ2_create3D(histogram_5_bins, number_of_bins, 1, 1, 32);
Ext.CLIJ2_histogram(image, histogram_5_bins, number_of_bins, 0, 5, false);
print("Histogram with " + number_of_bins + " bins");
Ext.CLIJ2_print(histogram_5_bins);

// make a histogram with 6 bins
number_of_bins = 6;
Ext.CLIJ2_create3D(histogram_6_bins, number_of_bins, 1, 1, 32);
Ext.CLIJ2_histogram(image, histogram_6_bins, number_of_bins, 0, 5, false);
print("Histogram with " + number_of_bins + " bins");
Ext.CLIJ2_print(histogram_6_bins);

image

Let me know if it helps! And furthermore, when you’re done with translating your whole workflow, let me know how fast it is using CLIJ2 compared to the original, if possible! :wink:

Thanks!

Cheers,
Robert