估计具有指数平滑和不规则事件的事件的发生率
abl*_*igh 14 algorithm poisson smoothing exponential
想象一下,我有一组x的测量值,这些测量值在时间t 0 ... t N被许多过程x 0 ... x N所取.假设在时间t我想基于我不知道的长期趋势的假设来估计x的当前值,并且可以从诸如指数平滑的算法预测x.由于我们有很多进程,N可以变得非常大,我不能存储多个值(例如以前的状态).
这里的一种方法是调整正常指数平滑算法.如果样品采取定期,我会保持一个估计Ÿ ñ这样的:
y n = α.y n-1 +(1 - α).X ñ
在采样不规则的情况下,这种方法并不是很好,因为许多样本在一起会产生不成比例的影响.因此,这个公式可以适用于:
ÿ Ñ = α Ñ.ý N-1 +(1 - α Ñ).X ñ
哪里
α Ñ = È -k(吨.ñ -吨N-1 )
IE根据前两个样本之间的间隔动态调整平滑常数.我很满意这种方法,似乎有效.这是这里给出的第一个答案,Eckner在2012年的论文(PDF)中给出了这些技术的一个很好的总结.
现在,我的问题如下.我想调整上述内容来估计发生率.偶尔会发生一个事件.使用类似的指数技术,我想估计事件发生的速率.
两个明显的策略是:
- 使用第一种或第二种技术,使用最后两个事件之间的延迟作为数据系列x n.
- 使用第一或第二技术,使用最后两个事件之间的延迟的倒数(即速率)作为数据序列x n.
据我所知,这些都不是好的策略.首先,采取每500毫秒(一方面)发生的事件和一个200毫秒延迟和另一方800毫秒延迟发生的事件.显然这些都是每秒发生两次,所以给出的速率估计应该是相同的.忽略上一个样本的时间似乎是蛮干的,所以我将专注于第二个策略.使用延迟(而不是倒数)并不是一个好的预测因子,因为模拟200ms/800ms样本流产生的估计值约为1.5(基于倒数的平均值不是平均值的倒数).
但是,更重要的是,这两种策略都不能应对实际发生的事情,即突然所有事件都会停止很长一段时间.因此,y的"最新"值是最后一个事件的值,因此永远不会计算出速率的估计值.因此,费率似乎是不变的.当然,如果我回顾性地分析数据,这不会是一个问题,但我正在实时分析它.
我意识到另一种方法是定期运行一些线程(例如每10秒)并计算这10秒间隔内的出现次数.由于统计数据不经常需要,因此我的资源非常繁重,而且我不愿意运行一个线程来调查因互斥问题导致的所有问题.因此,我想(以某种方式)使用一种算法,该算法通过(例如)自上一次采样以来的时间来调整读取的状态.这似乎是一种合理的方法,好像性能是在独立于样本选择的时间进行测量的,测量时间平均是样本之间的一半时间,所以非常粗略的未平滑估计的速率将是倒数的一半.自上次样本以来的时间.为了进一步复杂化,我的测量时间不会独立于样品.
我有一种感觉,这有一个简单的答案,但它是在逃避我.我有一种感觉,正确的路线是假设事件是泊松分布的,并根据自上次样本和某种形式的移动平均线以来的间隔得出λ的估计值,但是我的统计数据太生疏了,无法使其工作.
这里有一个近乎愚蠢的问题,但答案似乎并不十分令人满意(我希望我解释了原因).我补充说卡尔曼滤波器似乎是重量级的方法,因为我有一个变量来估计并且对它一无所知.还有其他一些近乎愚蠢的行为,其中大多数都建议保留较大的价值箱(从记忆的角度来看这里不现实)或者没有解决上述两个问题.
Ilm*_*nen 18
首先,如果您假设事件本身的发生率是恒定的(或者您只对其长期平均值感兴趣),那么您可以简单地将其估计为:
λ*= N /(t - t 0)
其中t是当前时间,t 0是观察的开始,N是自t 0以来观察到的事件的数量,并且λ*是真实频率λ的估计.
在这一点上,有用的是要注意上面给出的估算公式可以重新表述为积分:
λ*=积分(δ 事件(τ)dτ)/积分(1dτ)
其中积分τ的可变范围从吨0至吨,和δ 事件(τ)=总和(δ(τ - 吨我),我 = 1 .. Ñ)是总和Ñ的狄拉克δ函数,具有单个增量- 在每个事件的发生时间t i处的峰值i.
当然,这将是一种完全无用的计算 λ*的方法,但事实证明它是一个概念上有用的公式.基本上,查看此公式的方法是函数δ 事件(τ)测量事件数量在时间τ增加的瞬时速率,而第二个被积函数(仅为常数1)测量时间速率随着时间的推移而增加(当然,这只是每秒一秒).
好的,但是如果频率λ本身可能会随着时间的推移发生变化,并且您想要估计其当前值,或者至少估计其近期的平均值,该怎么办？
使用上面给出的积分比公式,我们可以简单地通过一些称重函数w(τ)对两个被积函数进行加权来获得这样的估计,该函数偏向近期:
λ*recent =积分(δ 事件(τ)w(τ)dτ)/积分(w(τ)dτ)
现在,剩下的就是选择合理的w(τ),使这些积分简化为易于计算的东西.事实证明,如果我们为某些衰减率k选择形式为w(τ)= exp(k(τ - t))的指数衰减加权函数,则积分可简化为:
λ*recent = sum(exp(k(t i - t)),i = 0 ... N)k /(1 - exp(k(t 0 - t)))
在极限为t 0 →-∞(即,实际上,当总观察时间(t - t 0)远大于权重衰减时间刻度1/k)时,这进一步简化为:
λ*recent = k sum(exp(k(t i - t)),i = 0 .. N)
唉,天真地应用这个公式仍然需要我们记住所有的事件时间t i.然而,我们可以使用与计算通常指数加权平均值相同的技巧 - 给定加权平均事件率λ*near(t')在某个较早时间t',并假设在t'和t之间没有发生新事件,我们可以简单地计算当前加权平均事件率λ*recent(t):
λ*recent(t)= exp(k(t' - t))λ*recent(t')
此外,如果我们现在观察到恰好在时间t发生的新事件,则事件之后的加权平均事件率变为:
λ*recent(t)= k + exp(k(t' - t))λ*recent(t')
因此,我们得到一个非常简单的规则:所有我们需要存储的时候牛逼最后先前观察到的事件,而最近估计事件发生率λ*持续刚才说的之后的事件.(我们可以初始化这些例如t last = t 0和λ*last = 0;实际上,在λ*last = 0的情况下,t last的值没有区别,尽管对于非零λ*最后它确实如此.)
每当发生新事件时(在时间t new),我们将这些值更新为:
λ*最后 ← ķ + EXP(ķ(吨持续 - Ť 新))λ*最后
吨最后 ← 吨新
每当我们想知道当前时间t的最近事件率时,我们只需将其计算为:
λ*(吨)= EXP(ķ(吨持续 - 吨))λ*最后
PS.为了校正对t last的(任意)初始值的初始偏差,我们可以加回1 /(1 - exp(k(t 0 - t)))校正项,当我们假设t时我们简化了» t 0.要做到这一点,只需从开始牛逼最后 = 0 牛逼 = 牛逼0,更新牛逼持续同上,但在时间计算估计,近期事件发生率平均牛逼为:
λ*更正件(吨)= EXP(ķ(吨持续 - 吨))λ*最后 /(1 - EXP(ķ(吨0 - 吨)))
(这里,t 0表示开始测量事件的时间,而不是第一个事件的发生时间.)
这将消除初始偏向零,代价是增加早期方差.这是一个显示校正效果的示例图,对于k = 0.1和真实平均事件率为2:
红线表示没有初始偏差校正的λ*(t)(从λ*(t 0)= 0开始),而绿线表示偏差校正估计λ*corr(t).
PPS.如上图所示,如上所计算的λ*不是时间的连续函数:每当事件发生时它会跳跃k,并且当事件不发生时它会指数衰减到零.
如果您更喜欢更平滑的估计,您可以计算λ*本身的指数衰减平均值:
λ**(t)=积分(λ*(τ)exp(k 2(τ - t))dτ)/积分(exp(k 2(τ - t))dτ)
其中λ*如上计算出的指数衰减平均事件速率,ķ 2是用于第二平均衰减率,并且积分超过-∞<τ≤ 吨.
此积分也可以通过上述逐步更新规则计算:
λ**最后 ← w ^(Δ 吨)λ*最后 + EXP( - ķ 2 Δ 吨)λ**最后
λ*最后 ← ķ 1 + EXP( - ķ 1 Δ 吨)λ*最后
吨最后 ← 吨新
其中ķ 1和ķ 2是用于在第一和第二平均值的衰变率,Δ 吨 = 吨新 - 吨最后是事件之间经过的时间,并且:
w ^(Δ 吨)= ķ 2(EXP( - ķ 2 Δ 吨) - EXP( - ķ 1 Δ 吨))/(ķ 1 - ķ 2)
如果k 1 ≠ k 2,或
w ^(Δ 吨)= ķ Δ 吨 EXP( - ķ Δ 吨)
如果ķ 1 = ķ 2 = ķ(来自前产生的作为时(极限后者表达ķ 1 - ķ 2)→0).
要计算任意时间点t的第二个平均值,请使用相同的公式:
λ**(吨)= w ^(Δ 吨)λ*最后 + EXP( - ķ 2 Δ 吨)λ**最后
除了Δt = t - t last.
如上所述,通过应用合适的时间相关比例因子,也可以对此估计进行偏差校正:
λ**corr(t)=λ**(t)/(1 - S(t - t 0))
哪里:
小号(Δ 吨)=(ķ 1 EXP( - ķ 2 Δ 吨) - ķ 2 EXP( - ķ 1 Δ 吨))/(ķ 1 - ķ 2)
如果k 1 ≠ k 2,或
小号(Δ 吨)=(1 + ķ Δ 吨)EXP( - ķ Δ 吨)
如果k 1 = k 2 = k.
下图显示了这种平滑效果.红色和绿色线如上所示显示λ*(t)和λ*corr(t),而黄色和蓝色线显示λ**(t)和λ**corr(t),如k 1 = 0.1 所计算(如上所述)和k 2 = 0.2:
