如何使用 SciPy Image 应用统一过滤器，其中不计算边界外的数据？答案

【问题标题】：How to apply a uniform filter using SciPy Image where the data outside the boundary is not tallied?如何使用 SciPy Image 应用统一过滤器，其中不计算边界外的数据？
【发布时间】：2019-05-18 06:33:03
【问题描述】：

我在这方面花费了很多时间，并且我知道如何通过对边界行/列进行切片和索引来手动执行此操作，但是 SciPy 必须有一种更简单的方法。

我需要将CVAL（在mode=constant 时填充边缘的值）设置为 NaN，但是，这将返回 NaN。

我会用代码和图来解释：

import numpy as np
from scipy import ndimage
m = np.reshape(np.arange(0,100),(10,10)).astype(np.float)

使用 SciPy ndimage uniform filter 使用 3x3 kernel 计算均值：

filter = ndimage.uniform_filter(m, size=3, mode='constant')
print(filter[1][1]) # equal to 11
print(filter[9][9]) # I need 93.5, however it gets 41.55 due to zeros

如您所见，第一个值为 11，这与预期的一样，但是，对于边界上的任何单元格，它将用零填充值（我也尝试了所有其他模式）。

这是我需要实现的（左） vs mode=constant 和 CVAL=0（默认 0）

【问题讨论】：

标签： python image-processing scipy ndimage

【解决方案1】：

一种简单的方法是使用Normalized Convolution：

import numpy as np
from scipy import ndimage
m = np.reshape(np.arange(0,100),(10,10)).astype(np.float)

filter = ndimage.uniform_filter(m, size=3, mode='constant')    # normal filter result

weights = ndimage.uniform_filter(np.ones(m.shape), size=3, mode='constant')
filter = filter / weights    # normalized convolution result

print(filter[1][1]) # equal to 11
print(filter[9][9]) # equal to 93.49999999999994 -- rounding error! :)

如果所有数据点都是 1 (weights)，我们会计算过滤器的结果。这显示了每个过滤器窗口中有多少数据元素，并在除边界附近以外的任何地方返回值 1，该值按比例减小。通过将过滤结果除以这些权重，我们将数据域之外的零考虑在内来校正平均。

【讨论】：

谢谢！我标记为正确，因为这非常简单且非常快速。

【解决方案2】：

这个建议不太理想，因为它比uniform_filter 慢，但它会做你想做的事。

使用nan作为常量值的想法，您可以使用ndimage.generic_filter而不是uniform_filter来实现统一过滤器，numpy.nanmean作为通用过滤器函数。

例如，这是您的示例数组m：

In [102]: import numpy as np

In [103]: m = np.reshape(np.arange(0,100),(10,10)).astype(np.float)

应用generic_filter，以numpy.nanmean为要应用的函数：

In [104]: from scipy.ndimage import generic_filter

In [105]: generic_filter(m, np.nanmean, mode='constant', cval=np.nan, size=3)
Out[105]: 
array([[ 5.5,  6. ,  7. ,  8. ,  9. , 10. , 11. , 12. , 13. , 13.5],
       [10.5, 11. , 12. , 13. , 14. , 15. , 16. , 17. , 18. , 18.5],
       [20.5, 21. , 22. , 23. , 24. , 25. , 26. , 27. , 28. , 28.5],
       [30.5, 31. , 32. , 33. , 34. , 35. , 36. , 37. , 38. , 38.5],
       [40.5, 41. , 42. , 43. , 44. , 45. , 46. , 47. , 48. , 48.5],
       [50.5, 51. , 52. , 53. , 54. , 55. , 56. , 57. , 58. , 58.5],
       [60.5, 61. , 62. , 63. , 64. , 65. , 66. , 67. , 68. , 68.5],
       [70.5, 71. , 72. , 73. , 74. , 75. , 76. , 77. , 78. , 78.5],
       [80.5, 81. , 82. , 83. , 84. , 85. , 86. , 87. , 88. , 88.5],
       [85.5, 86. , 87. , 88. , 89. , 90. , 91. , 92. , 93. , 93.5]])

【讨论】：

啊，是的通用过滤器，谢谢！我会检查性能。通常我会处理形状为 (1e6, 1e6) 的数组。
我怀疑你不会高兴的。 ?
是的，你是对的，这非常慢。 m = np.random.rand(1000,1000 统一过滤器约为 0.01 秒，而通用过滤器在我的计算机上约为 22 秒。
可能最好的办法是在完整图像上使用uniform_filter，然后再添加几行来修复边界。从您在问题中写的内容来看，我认为您已经在这样做了。