贝叶斯统计:观测者偏移

最新推荐文章于 2021-05-31 14:46:14 发布

duskwaitor

最新推荐文章于 2021-05-31 14:46:14 发布

阅读量1.1k

点赞数

本文链接：https://blog.csdn.net/duskwaitor/article/details/41549579

版权

文章作者说，由于工作关系，他要乘坐地铁，他们那里大约8分钟一趟地铁，如果他到地铁站时候发现人很少，他估计车刚走，得等7分钟左右，如果人很多，他认为很快就来车了，如果人多的不行了，估计交通阻塞，那么他就会打的。

然后作者就想到了，能否用贝叶斯估计来帮忙预测等车时间，然后再帮他决定是继续等地铁还是打的。

首先，做一个假设，旅客到达车站的时间分布是泊松过程，以lambda个人每分钟的节奏到达，lambda是个常量。

另外，对于地铁来说它的到达过程并非泊松过程，从Alewife出发到达波士顿的地铁，在高峰期每7-8分钟一趟，但当他们到达Kendall广场时，相邻列车的时间相隔是3-12分钟。于是我写了个脚本(我指的是原作者)，抽取15天中，到达kendall广场的南向列车的间隔时间分布图，如图1所示：

图1 列车时间间隔分布情况（紫色是真实分布，浅蓝是观察到的分布）

当你站在早晨4点到6点的月台上观察列车间隔时间时，你将看到这个分布，然而，如果你随机到达月台，你将看到不一样的分布，列车间隔的平均时间，要远高于真实的平均时间。

为毛？因为乘客更可能在时间大的区间到达车站，举一个简单的例子，5分钟一趟和10分钟一趟的车等可能到达，那么平均到达时间是7.5分钟，但10分钟间隔内到达车站的旅客，其数量就是5分钟间隔的两倍，如果我们调查，就会发现2/3的乘客10分钟间隔到达，1/3的乘客5分钟间隔到达，所以在一个乘客看来，地铁平均间隔时间是8.33分。

这就是观察者偏移，这种情况发生的机会很多，比如乘客总认为飞机装载的比实际满，班级的学生总认为班级人数多。

因而给定一个实际的地铁间隔分布，我们可以计算出乘客所看到的地铁间隔分布，BiasPmf就是做这类计算的：

def BiasPmf(pmf):
     new_pmf = pmf.Copy()
     for x,p in pmf.Items()
          new_pmf.Mult(t,x)
     new_pmf.Normalize()
     return new_pmf

等待时间

等待时间，y，就是顾客到达时间和下次列车到达时间的间隔，拖延时间，x，就是上次列车到达时间和旅客到达时间的间隔，zb = x+y.

给出zb的分布，我们就能够得出y的分布，从一个小例子讲起：

还是这个例子，2/3的乘客10分钟间隔到达，1/3的乘客5分钟间隔到，而y在5分钟和10分钟的间隔内均匀分布，所以整个分布就是两个均匀分布的混合分布，每个均匀分布的权重正比于每个间隔出现的概率。

下面代码就根据zb的分布，来计算y的分布：

def PmfOfWaitTime(pmf_zb):
      metapmf = thinkbayes.Pmf()
      for gap,prob in pmf_zb.Items()
           uniform = MakeUniformPmf(0,gap)
           metapmf.Set(uniform,prob)
      pmf_y = thinkbayes.MakeMixture(metapmf)
      return pmf_y

现在我们得到z,zb,y的cdf图形如下：

图2 z,zb,y的cdf图(均值：z 7.8,zb 8.8 y 4.4)

预测等待时间

好了，让我们接下来来预测等待时间吧。

实际情景是这个样子，有一天，我到了地铁站，发现站台上有10个人，那么我需要等多长时间才能等到下一班车？按照往常的套路，让我们先从简单入手，假设我们事先已经知道了z分布，而且我们也知道了lambda，每分钟两个人进入地铁站，在这种情况下：

1.使用z的分布计算zp的先验分布，zp=x+y，也就是乘客看到的列车间隔时间

2.我们使用候车乘客人数估计x的分布，拖延时间

3.最后，我们使用y=zp-x来获取y的分布

还是介绍几个函数，具体意思不表

class Elapsed(thinkbayes.Suite):
    """Represents the distribution of elapsed time (x)."""

    def Likelihood(self, data, hypo):
        """Computes the likelihood of the data under the hypothesis.

        Evaluates the Poisson PMF for lambda and k.

        hypo: elapsed time since the last train
        data: tuple of arrival rate and number of passengers
        """
        x = hypo
        lam, k = data
        like = thinkbayes.EvalPoissonPmf(k, lam * x)
        return like


def PredictWaitTime(pmf_zb, pmf_x):
    """Computes the distribution of wait times.

    Enumerate all pairs of zb from pmf_zb and x from pmf_x,
    and accumulate the distribution of y = z - x.

    pmf_zb: distribution of gaps seen by random observer
    pmf_x: distribution of elapsed time
    """
    pmf_y = pmf_zb - pmf_x
    pmf_y.name = 'pred y'
    RemoveNegatives(pmf_y)
    return pmf_y

def EvalPoissonPmf(k, lam):
    """Computes the Poisson PMF.

    k: number of events
    lam: parameter lambda in events per unit time

    returns: float probability
    """
    # don't use the scipy function.  for lam=0 it returns NaN;
    # should be 0.0
    return scipy.stats.poisson.pmf(k, lam)

估计乘客到达率