The currently-used method of adaptive thresholding (binning the image in chunks) results in object artifacts. A better method would be to take the binned thresholds and perform a 2D cubic interpolation. It seems to be working in my implementation, which places the block threshold in the center of the bin.

Add this code to the end of the get_adaptive_threshold subroutine in threshold.py (the first section is identical to what already exists except that I am storing the block thresholds instead of adding them to thresh_out:

[code]
import scipy as sp
…]
#
# Loop once per block, computing the “global” threshold within the
# block.
#
block_threshold = np.zeros([nblocks[0],nblocks[1]])
for i in range(nblocks[0]):
i0 = int(i*increment[0])
i1 = int((i+1)increment[0])
for j in range(nblocks[1]):
j0 = int(j
increment[1])
j1 = int((j+1)increment[1])
block = image[i0:i1,j0:j1]
block_threshold
= get_global_threshold(threshold_method,
object_fraction = object_fraction,
two_class_otsu = two_class_otsu,
use_weighted_variance = use_weighted_variance,
assign_middle_to_foreground = assign_middle_to_foreground)

``````#
# Use a cubic spline to blend the thresholds across the image to avoid image artifacts
#
xStart = int(increment[0] / 2)
xEnd = int((nblocks[0] + 0.5) * increment[0])
yStart = int(increment[1] / 2)
yEnd = int((nblocks[1] + 0.5) * increment[1])

block_x_coords = np.linspace(xStart,xEnd, nblocks[0])
block_y_coords = np.linspace(yStart,yEnd, nblocks[1])
thresh_out_x_coords = np.linspace(0, int((nblocks[0] + 1) * increment[0]), thresh_out.shape[0])
thresh_out_y_coords = np.linspace(0, int((nblocks[1] + 1) * increment[1]), thresh_out.shape[1])

return thresh_out
``````

[/code]*

Hi,
Thanks for drafting this code. This issue was actually discussed briefly on the cp-dev mailing list, and we reached the same conclusion re: implementation. We have added a note to include this improvement for the next release.
-Mark

Thanks for posting this. I’ve pretty much taken your patch (svn rev 11623), just a few subtle modifications:

The spine order needs to be less than the number of blocks in order to satisfy an assertion in the interpolation package and generally will be. If not, I use a quartic or linear interpolation.

I noticed that the spline function has some artifacts for the blocks around the edge - it needs to be told that the bounding box for the spline function is from 0 to the end of the image, otherwise it uses a bounding box from the midpoint of one end block to the midpoint of the other.

I made a couple mods related to pixel positioning which are probably pedantic: the linear space starts at .5 and ends at the image size - .5 so that each member of the destination linear space is one pixel wide.

Hope you don’t mind me butting in with my changes - if they don’t work correctly, by all means get back in touch with me.

–Lee

I’ve been away from working with CellProfiler for a while so I didn’t check back to see your reply. I don’t mind the butting-in at all! I’m fairly novice at this compared to you guys, and your mods make perfect sense to me.