python - 使用 scipy generic_filter 和 numpy Median_filter 计算移动中位数会给出不同的输出

Question

我希望实现一个快速移动的中位数，因为我必须为我的程序做很多中位数。我想使用 python 内置函数，因为它们比我能做的更优化。

我的中位数应该做:

• 提取 5 个值，
• 删除中间的一个，
• 求剩余 4 个值的中位数。

基本上多次调用:

numpy.median(np.array([0, 1, 2, 3, 4])[np.array([True, True, False, True, True])])
# (1. + 3.) / 2. = 2.0

我发现了两个函数:scipy generic_filter 和 scipy Median_filter。我的问题是 generic_filter 提供了正确的输出，而不是median_filter，即使它们似乎具有相同的参数。而且，generic_filter比median_filter慢。所以我想知道我在调用median_filter时做错了什么，并使用这个来提高速度。

import numpy as np
import scipy.ndimage as sc

v = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10])

print(sc.generic_filter(v, sc.median, footprint=np.array([1, 1, 0, 1, 1]), mode = "mirror", output=np.float64))
%timeit sc.generic_filter(v, sc.median, footprint=np.array([1, 1, 0, 1, 1]), mode = "mirror", output=np.float64)

print(sc.median_filter(v, footprint=np.array([1, 1, 0, 1, 1]), output=np.float64, mode="mirror"))
%timeit sc.median_filter(v, footprint=np.array([1, 1, 0, 1, 1]), output=np.float64, mode="mirror")

如您所见，generic_filter 给出了正确的输出:[1.5 1.5 2.3.4.5.6.7.8.8.5 8.5]每个循环 327 µs ± 15.2 µs(7 次运行的平均值 ± 标准偏差，每次 1000 个循环)

和median_filter更快，但我不明白它的输出:[2. 2.3.4.5.6.7.8.9.9.9.]每个循环 12.4 µs ± 217 ns(7 次运行的平均值 ± 标准差，每次 100000 个循环)

你知道我的电话有什么问题吗？

Answer 1

唯一的区别似乎是由于如何处理“关系”:

• sc.median 返回关系的平均值
• sc.median_filter 似乎系统地返回较大的值

给定方式median_filter is implemented对于“偶数个元素的中位数应该返回关系的平均值”的情况，有效处理特殊/特定是很尴尬的

我已经编写了一个可以处理这种情况的版本:

from scipy.ndimage.filters import _rank_filter

def median_filter(input, footprint, output=None, mode="reflect", cval=0.0, origin=0):
    filter_size = np.where(footprint, 1, 0).sum()
    rank = filter_size // 2
    result = _rank_filter(
        input, rank, None, footprint, output, mode, cval, origin, 'dummy')
    if filter_size % 2 == 0:
        if result is output:
            tmp = result.copy()
        else:
            tmp = result
        rank -= 1
        assert rank > 0
        result = _rank_filter(
            input, rank, None, footprint, output, mode, cval, origin, 'dummy')
        # fix up ties without creating any more garbage
        result += tmp
        result /= 2
    return result

但它有点笨重，并且使用 scipy 的内部功能(我使用的是 1.3.0)，因此将来可能会崩溃

在我的机器上，这些基准测试为:

• sc.generic_filter 每个循环需要 578 µs ± 8.51 µs
• sc.median_filter 每个循环需要 27.4 µs ± 1.37 µs
• 我的 median_filter 每次循环需要 65.6 µs ± 1.29 µs