Code audit: usage of `img_as_float`

[jni]

There's been lots of discussion about our automatic range conversions, especially when it comes to img_as_float. (See for example #3009 and #3052.) In #3052 @stefanv suggested that we do an audit to check which functions will be affected if we don't ensure that floats are in the range [0, 1]. I went through every non-testing use of img_as_float in our code base to check. Below, I've divided the functions into three cases: (1) function will not work correctly if image is out of range; (2) function will return a correspondingly-scaled image (ie it will work correctly, but the output is not scale invariant), and (3) functions where the output will be unchanged (ie we can safely remove range conversion altogether).

Range essential for correctness

• exposure.intensity_range will return [0, 1] for range_values='dtype' or range_values='float'. However this behaviour in no way depends on img_as_float, let alone on img_as_float raising an error.
• feature.draw_haar_like_features does need the image to be in the correct range because colours are normalised to be in [0, 1]
• Same for feature.draw_multiblock_lbp.
• feature.plot_matches needs the range to make sure that the random color line drawn over it makes sense.
• restoration.denoise_wavelet probably needs a specific range, as it converts images to ycbcr and color conversion would assume a certain range for RGB input.
• segmentation.felzenszwalb method divides scale by 255 to be consistent with the reference implementation, I presume by assuming that the image values are in [0, 1] instead of [0, 255]. Therefore, here we would need to be careful to make sure that the scale matches the input image scale.
• segmentation.mark_boundaries "paints" colour onto the image, so the input image needs to be in [0, 1]. (Note that it currently doesn't check whether a float image is in that range, it just uses img_as_float; so this needs fixing, as do other cases of drawing covered in this list.)
• future.rag.show_rag needs the [0, 1] magnitude to display correctly in mpl

Corresponding magnitude change in output

• feature.corner.hessian_matrix_det and hessian_matrix_elems use it. There will be a difference in magnitude if the image is out of range but everything should still "just work".
• Same for feature.corner_moravec
• feature.structure_tensor as above.
• feature.CENSURE -- not sure about this one, but it has some thresholding which probably needs to change if we change the scale of the image. Uses _prepare_grayscale_input_2D, a convenience function that calls img_as_float.
• feature.corner_fast uses a threshold to compare to image sums, so threshold choice is dependent on image range. Uses _prepare_...
• feature.ORB uses corner_fast and brief. Brief is not sensitive to image range (see below) but fast corners are.
• filters.edges.{sobel,prewitt,scharr} and variants
• feature.blob_doh default threshold is intended for floats scaled in range [0, 1]
• restoration.denoise_bilateral takes a parameter sigma that depends on the range of the image, which currently is converted with img_as_float. So the sigma parameter would need to be scaled depending on the scale of the image values.
• same for denoise_total_variation; however in this case, the default is to estimate the sigma from the data, so it should work normally out of the box.
• quickshift and SLIC cluster on colour and coordinate space, so the relative magnitude parameter needs to change if the scale of the image range changes. The default values for these are rarely what you want, though, so functionally these functions should be mostly unaffected by the image range. (Though of course existing pipelines would break if we stopped range-converting.)
• Same for active_contours.
• I wasn't sure about random_walker... @JDWarner?

Behaviour unchanged / totally insensitive to image magnitude

• feature.daisy uses img_as_float but does not depend in any way on the range of the image, because the detector compares image magnitudes at different angles, so only the relative values of pixels to other pixels matters. (Incidentally this code needs a good cleanup!)
• feature.blob_* use img_as_float followed by some of the scale-sensitive functions mentioned above, but they then use local maxima to find peaks, so they don't depend on the magnitude at all. (one exception is the threshold parameter for blob_doh , for which an appropriate default is dependent on the magnitude of the floating point data).
• feature.brief only checks for relative magnitude between keypoints so behaviour is unchanged by range.
• feature.corner_orientation uses atan so is scale invariant
• filters.gaussian will not be affected (already lets values through when a float image is passed; otherwise converts when preserve_range is False (default)).
• restoration.inpaint_biharmonic (I think). Uses sparse linear algebra to solve an equation based on the data. No assumed data values that I can tell.
• transform.warps (values will be passed through as expected.)
• transform.pyramids (as above)
• transform.seam_carving (only cares about relative magnitudes)

Other notes

• def extract in feature/brief contains both assert_nD and _prepare_...2D, which is redundant.
• same for def detect in feature/censure.
• "minimas" and "maximas" are not correct ("minima" and "maxima" are already plural) in feature.censure.
• filters.unsharp_mask has some magic to detect negative values. img_as_float is used, which does not do range checking for float images, so this could result in buggy behaviour. (Although the preserve_range keyword mitigates that.)

The above comes from me examining the code. A more robust audit would use @emmanuelle's API crawl to test images on each function, then test 2x image and check how the output changes:

• not at all
• constant factor
• more complex change

Of course, we would still miss any functions requiring more than one input.

[scottstanie]

Is there any plan for an intermediate fix to allow the functions under "Behaviour unchanged / totally insensitive to image magnitude" to use float32 unscaled images?

for instance with feature.blob_*, right now it's necessary to do something like image.astype('float64'), which makes the img_as_float as no-op and avoids the ValueError.

[jni]

@scottstanie if I'm not mistaken, #3052 would allow this, is this correct? I'm hoping to get that merged later this week... :crossed_fingers:

[scottstanie]

Yep! hadn't seen that but it would fix it exactly as I'd hoped (making the float32 same as float64)

[emmanuelle]

@jni random walker was written to be insensitive to image magnitude, the beta parameter is normalized by img.std() but I just saw that this is a bug, beta should be normalized by the variance and not the standard deviation. Anyway, I think random walker falls in the category of algorithms (like SLIC) where you have to tune a parameter to get what you want.

[grlee77]

Summarizing from above

Below, I used rescale_to_float to refer to the function currently named img_as_float (see #6318).

Parameter Scaling

For the majority of functions currently using rescale_to_float, removing this just results in an output range that is different by exactly the former rescaling factor. For example, gaussian filtering a uint8 image without rescaling by 1/255.0 will result in an output that is 255 times larger than previously. The functions in skimage.transform and skimage.filters tend to fall in this category.

There are a few exceptions to this behavior to consider:

• Functions where the output is a binary segmentation or label image. In this case the output range is invariant to image scale although some keyword parameters to functions may need adjustment as detailed in the "Parameter Scaling" table below.

• There are a few functions where the output is proportional to the square of the amplitude of the input. These include structure_tensor, hessian_matrix and hessian_matrix_det. This seems okay.

• A handful of functions clip output to range [0, 1]. This needs to either be removed or the rescaling should be preserved. For example, color.combine_stains clips the output to [0, 1]

• exposure.adjust_gamma, exposure.adjust_log, exposure.adjust_sigmoid and exposure.is_low_contrast rely on use of dtype_limits. In this case, it seems that floating point inputs need to already be in range [0, 1], but this is not currently being checked or enforced.

• Color conversions need to be audited individually. A number of these assume normalization to floats in [0, 1], but a subset (e.g. rgb2gray) are safe to apply without any normalization.

Visualization functions are one case where we should continue to normalize to [0, 1]. These include feature.plot_matches, feature.draw_haar_like_features, feature.draw_multiblock_lbp, future.rag.show_rag, segmentation.mark_boundaries, that involve visualization with Matplotlib's imshow. Most or all of these currently rescale int inputs, but not float ones, so that should be fixed!

Parameter Scaling

The following table lists functions that have internal use of rescale_to_float, but where appropriate choice of one of the keyword-arguments depends on the absolute image scale. The table below shows the relative scaling that needs to be applied to a parameter when an image is scaled by image = a * image.

function	parameter name(s)	scaling
feature.blob_log	threshold	a
feature.blob_dog	threshold	a
feature.blob_doh	threshold	a
feature.corner_fast	threshold	a
feature.CENSURE	non_max_threshold	a
feature.ORB	fast_threshold	a
feature.SIFT	c_dog?, c_edge?	a
restoration.denoise_tv_bregman	weight	1 / a
restoration.denoise_tv_chambolle	weight	a
restoration.denoise_bilateral	sigma_color	a
restoration.denoise_wavelet	sigma (but default sigma=None estimation is scale insensitive)	a
segmentation.felzenszwalb	sigma	a
segmentation.random_walker	beta, tol?	a?
util.random_noise	var	a**2

Notes:

• feature.blob_* also have a threshold_rel that can be used instead and is insensitive to absolute scale.

• Prior to v0.19, for SLIC superpixels (segmentation.slic), the choice of compactness was image amplitude dependent. As of v0.19, internal normalization to range [0, 1] is always performed, even for floating point dtypes.

• Apparent bug in the segmentation.felzenswalb Cython code. It calls rescale_to_float, but then separately uses scale /= 255. Thus, the scale is reduced even in cases such as floating point input that aren't rescaled!

• For restoration.denoise_wavelet in addition to an internal rescale_to_float for integer images. Clipping of the output to [0, 1] or [-1, 1] is automatically applied only when the input is of integer type.

[alexdesiqueira]

Thank you for #6318 , @grlee77; it looks great already! I was thinking, though... wouldn't `rescaletogive a different idea on the image context? Couldn't we useconvertto*` instead...? :slightly_smiling_face:

Edit: I asked the same in the past. Any ideas?