我在寻找一种快速有效的方法来计算一组数据的稳健的移动尺度估计。我使用的是1d数组,通常是3-400k个元素。直到最近,我一直在处理模拟数据(没有灾难性的异常值),优秀的瓶颈包中的move_std函数对我很有用。然而,由于我已经过渡到嘈杂的数据,std不再表现得足够好用了。在
在过去,我使用了一个非常简单的双权中方差代码元素来处理不良分布的问题:
def bwmv(data_array):
cent = np.median(data_array)
MAD = np.median(np.abs(data_array-cent))
u = (data_array-cent) / 9. / MAD
uu = u*u
I = np.asarray((uu <= 1.), dtype=int)
return np.sqrt(len(data_array) * np.sum((data_array-cent)**2 * (1.-uu)**4 * I)\
/(np.sum((1.-uu) * (1.-5*uu) * I)**2))
不过,我现在使用的数组足够大,速度非常慢。有没有人知道提供这种估计器的软件包,或者对如何快速有效地处理这个问题有什么建议?在
我在类似的情况下使用了一个简单的低通滤波器。在
从概念上讲,您可以使用
fac = 0.99; filtered[k] = fac*filtered[k-1] + (1-fac)*data[k]
获得平均值的移动估计值,这是非常有效的实现(在C中)。一个比这个稍微更花哨的IIR滤波器,butterworth低通,很容易在scipy中设置:要获得“量表”的估计值,可以从数据中减去这个“平均估计值”。这实际上将低通滤波器转换为高通滤波器。把它的abs()放到另一个低通滤波器中。在
结果可能如下:
完整脚本:
^{pr2}$显然,butter()参数需要针对您的问题进行调整。如果您将order设置为1而不是2,那么您将得到我首先描述的简单过滤器。在
免责声明:这是工程师对问题的看法。这种方法在统计学或数学上可能都不合理。另外,我不确定它是否真的解决了你的问题(如果解决不了,请解释得更好),但别担心,我做这件事很有意思,不管怎样;-)
相关问题 更多 >
编程相关推荐