0

EDIT: Problem solved in comments, thank you for your help.
EDIT 2: I think that the linked solution is mathematically not correct, but maybe I am wrong. The problem is that you divide by the window size but you do not add the quadratic difference together. Solution: My solution now is to calculate the precision by np.finfo(np.float).precision and set every value to zero, that is below pow(1, np.finfo(np.float).precision)] because that should be the upper bound for the machine epsilon error. If you have negative values in your matrix, you have to multiply them by -1 before you check if the value is below your error threshold. Otherwise values like -2 would be set to 0 as well.

I have an 2D array with intensity values derived from an image (52x111). To calculate the local standard deviation I used the answer from another thread (improving code efficiency: standard deviation on sliding windows)
The smallest intensity value is 0.0, as parameters for the filter function I used mode = "constant" and cval = 0.

My mean filter looks like this:

mean = uniform_filter(img_data, size=3, cval= 0, mode="constant", origin=0)

The result of np.amin(mean) is -2.47949808833e-15.

How is it possible, that the filter yield to negative values?

To clarify things, the whole relevant code looks like this:

mean = uniform_filter(img_data, size=3, cval= 0, mode="constant", origin=0) 
mean_of_squared = uniform_filter(img_data**2, 3, cval= 0, mode='constant', origin=0) 
squared_mean = mean * mean 
stdev = (mean_of_squared - squared_mean) ** 0.5

As alternative I use the concolve function from scipy

window_size = 3
mean_filter = np.ones((window_size, window_size)) / 9
mean = convolve(img_data, mean_filter)

which works fine, but is much slower than the uniform_filter approach.

As far as I remember all of the coordinates that I inspected manually are just windows of 0.0, but I am not sure if this is the general case. Any ideas why uniform_filter is acting like it does?

Community
  • 1
  • 1
Xenthor
  • 423
  • 6
  • 12
  • Possible duplicate of [Looking for a specific task similar to sci.py's "ndimage.filters.uniform\_filter"](http://stackoverflow.com/questions/37835238/looking-for-a-specific-task-similar-to-sci-pys-ndimage-filters-uniform-filter) – Benjamin Jun 24 '16 at 13:43
  • @Benjamin I think you got me wrong. This is not a duplicate. What you posted is what I try to do. But the sqrt operation is not possible because the filter yield to negative values. With your formular this could be solved, because you square everything, but it still does not make any sense, that the filter yield to negative values if all my original values are above 0. – Xenthor Jun 24 '16 at 14:04
  • Can you clarify "the filter yield to negative values"? Taking the square difference in the SD should ensure that you are not taking the square root of a negative number. – Benjamin Jun 24 '16 at 14:07
  • I think I add more code to clarify it. I do the following: `mean = uniform_filter(img_data, size=3, cval= 0, mode="constant", origin=0)` `mean_of_squared = uniform_filter(img_data**2, 3, cval= 0, mode='constant', origin=0)` `squared_mean = mean * mean` `stdev = (mean_of_squared - squared_mean) ** 0.5` Mathematically this should be fine and every value is squared at some point in this calculation, but i still get values like `-2.47949808833e-15` in my stdev matrix. – Xenthor Jun 24 '16 at 14:08
  • @Xenthor Please add that code to the question. It is probably important. – Warren Weckesser Jun 24 '16 at 14:32
  • Not clear, but it looks like you are assuming `(a-b)**2 == a**2 - b**2`, which is not the case. The solution I list in the duplicate is correct. – Benjamin Jun 24 '16 at 14:32
  • No I assume that Var(X) = E(X^2) - E(X)^2 and this is correct as well. Sorry but I feel not comfortable just using your solution, while using a function which gives me suspicious values, but should give me correct ones from the mathematical point of view. – Xenthor Jun 24 '16 at 14:47
  • 3
    In that case, your very small value is probably a zero, very close to `np.finfo(np.float).eps`. It's not an issue with the uniform filter... just set them to zero before taking the sqrt, or use a mask to apply it only to non-zero values (or use my approach which will avoid the issue all together). See wikipedia article for variance: "On computational floating point arithmetic, this equation should not be used, because it suffers from catastrophic cancellation if the two components of the equation are similar in magnitude." – Benjamin Jun 24 '16 at 15:10
  • @Benjamin Thank you very much, I think this is very plausible and explains everything. I nearly forgot, that this problem can arise. – Xenthor Jun 25 '16 at 06:59

0 Answers0