Using UniformDiscrete and DensityDist for Mixture distribution throws IndexError: axis 1 is out of bounds [-1, 1)
Using UniformDiscrete and DensityDist for Mixture distribution throws IndexError: axis 1 is out of bounds [-1, 1)
我有以下模型，其中使用 Jeffrey 的先验几何分布作为其中一种分布。我有一个派生自其他分布的分布，并在混合分布中使用它。我得到
IndexError: axis 1 is out of bounds [-1, 1) error.
请在下面找到我的代码。
actual_visit_mod_90_days = df_visit_period['mod_90_days'].values
def jeffreys_logp(theta):
#log(((1-theta)**-1) * theta**-0.5 )
return -(T.log(1-theta) + 0.5 * T.log(theta))
with Model() as pattern_study:
dist1 = Binomial('dist1',n=14,p=0.5)
dist2 = DiscreteUniform('dist2',lower = 0, upper = 89)
pp = pm.DensityDist('pp',jeffreys_logp, testval = 0.05)
dist3 = Geometric('dist3',pp)
someDervDist = pm.DensityDist('someDervDist', lambda x: dist2-dist1+dist3 -1)
likelihoodDist1 = DiscreteUniform.dist(lower=0, upper=89)
likelihoodDist2 = pm.DensityDist.dist(lambda x: pm.logsumexp(someDervDist))
p=0.05
p_vd = Mixture('p_vd',[p,1-p],[likelihoodDist1,likelihoodDist2],observed = actual_visit_mod_90_days)
错误信息是
---------------------------------------------------------------------------
IndexError Traceback (most recent call last)
<ipython-input-19-1c9d266976ce> in <module>()
28 p=0.05
29
---> 30 p_vd = Mixture('p_vd',[p,1-p],[likelihoodDist1,likelihoodDist2],observed = actual_visit_mod_90_days)
31
32 #print pm.distributions.generate_samples(postponeProbability)
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\distributions\distribution.pyc in __new__(cls, name, *args, **kwargs)
35 total_size = kwargs.pop('total_size', None)
36 dist = cls.dist(*args, **kwargs)
---> 37 return model.Var(name, dist, data, total_size)
38 else:
39 raise TypeError("Name needs to be a string but got: {}".format(name))
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\model.pyc in Var(self, name, dist, data, total_size)
780 var = ObservedRV(name=name, data=data,
781 distribution=dist,
--> 782 total_size=total_size, model=self)
783 self.observed_RVs.append(var)
784 if var.missing_values:
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\model.pyc in __init__(self, type, owner, index, name, data, distribution, total_size, model)
1228
1229 self.missing_values = data.missing_values
-> 1230 self.logp_elemwiset = distribution.logp(data)
1231 # The logp might need scaling in minibatches.
1232 # This is done in `Factor`.
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\distributions\mixture.pyc in logp(self, value)
111 w = self.w
112
--> 113 return bound(logsumexp(tt.log(w) + self._comp_logp(value), axis=-1).sum(),
114 w >= 0, w <= 1, tt.allclose(w.sum(axis=-1), 1),
115 broadcast_conditions=False)
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\distributions\mixture.pyc in _comp_logp(self, value)
83 except AttributeError:
84 return tt.stack([comp_dist.logp(value) for comp_dist in comp_dists],
---> 85 axis=1)
86
87 def _comp_means(self):
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\theano\tensor\basic.pyc in stack(*tensors, **kwargs)
4583 dtype = scal.upcast(*[i.dtype for i in tensors])
4584 return theano.tensor.opt.MakeVector(dtype)(*tensors)
-> 4585 return join(axis, *[shape_padaxis(t, axis) for t in tensors])
4586
4587
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\theano\tensor\basic.pyc in shape_padaxis(t, axis)
4475 if not -ndim <= axis < ndim:
4476 msg = 'axis {0} is out of bounds [-{1}, {1})'.format(axis, ndim)
-> 4477 raise IndexError(msg)
4478 if axis < 0:
4479 axis += ndim
IndexError: axis 1 is out of bounds [-1, 1)
对我来说，问题似乎出在 DensityDist.dist 函数上。我也不确定我是否遗漏了什么。
编辑更新：考虑到事件 1 的历史模式 (actual_visit_mod_90_days)，我正在尝试了解事件 1 从开始日期起多少天后会发生。最初，我将 someDervDist 定义为一个等式：
someDervDist = dist2-dist1+dist3-1
给出了决定事件 1 的事件的各种概率之间的关系。
likelihoodDist2 = pm.DensityDist.dist(lambda x: pm.logsumexp(someDervDist))
会给出 logsumexp(someDervDist.logp()) 我没有写成 .logp 因为我认为它会给出对数概率（嗯，我在这里错了吗？）但是，我意识到这需要日志概率并将 someDervDist eqn 更改为 DensityDist 分布（我也尝试使用 Deterministic）。
然后，我假设这偏离事件 1 的周期性发生模型的概率很小 p=0.05。基本上我正在寻找类似 Mixture([p, 1-p], [ DiscreteUniform(lower = 0, upper = 89), someDervDist]).
更新的答案，毕竟使用 Constant 发行版不是一个好主意，它很快就会被弃用。所以也许一个小技巧会有所帮助，只需使用均值 someDervDist 和小方差的分布。
actual_visit_mod_90_days = [29, 76, 66, 66, 10, 63, 69, 75, 66, 87, 58]
with pm.Model() as pattern_study:
dist1 = pm.Binomial('dist1', n=14, p=0.5)
dist2 = pm.DiscreteUniform('dist2', lower = 0, upper = 89)
pp = pm.DensityDist('pp', jeffreys_logp, testval = 0.5)
dist3 = pm.Geometric('dist3', pp)
someDervDist = pm.Deterministic('someDervDist', dist2-dist1+dist3-1)
likelihoodDist2 = pm.DiscreteUniform.dist(someDervDist-1, someDervDist + 1)
likelihoodDist1 = pm.DiscreteUniform.dist(lower=0,upper=89)
p= 0.05
p_vd = pm.Mixture('p_vd', [p, 1-p],
[likelihoodDist1, likelihoodDist2],
observed = actual_visit_mod_90_days)
trace = pm.sample(1000, step=pm.Metropolis())
作为替代方案，您可以尝试使用模型的连续版本，如下所示：
with pm.Model() as pattern_study:
dist1 = pm.Normal('dist1', 7, 1.5)
dist2 = pm.Uniform('dist2', lower = 0, upper = 89)
dist3 = pm.Exponential('dist3', 1)
someDervDist = pm.Deterministic('someDervDist', dist2-dist1+dist3-1)
likelihoodDist2 = pm.Normal.dist(someDervDist, 0.1)
likelihoodDist1 = pm.Uniform.dist(lower=0,upper=89)
p = 0.05
p_vd = pm.Mixture('p_vd', [p, 1-p],
[likelihoodDist1, likelihoodDist2],
observed = actual_visit_mod_90_days)
trace = pm.sample(1000)
我有以下模型，其中使用 Jeffrey 的先验几何分布作为其中一种分布。我有一个派生自其他分布的分布，并在混合分布中使用它。我得到
IndexError: axis 1 is out of bounds [-1, 1) error.
请在下面找到我的代码。
actual_visit_mod_90_days = df_visit_period['mod_90_days'].values
def jeffreys_logp(theta):
#log(((1-theta)**-1) * theta**-0.5 )
return -(T.log(1-theta) + 0.5 * T.log(theta))
with Model() as pattern_study:
dist1 = Binomial('dist1',n=14,p=0.5)
dist2 = DiscreteUniform('dist2',lower = 0, upper = 89)
pp = pm.DensityDist('pp',jeffreys_logp, testval = 0.05)
dist3 = Geometric('dist3',pp)
someDervDist = pm.DensityDist('someDervDist', lambda x: dist2-dist1+dist3 -1)
likelihoodDist1 = DiscreteUniform.dist(lower=0, upper=89)
likelihoodDist2 = pm.DensityDist.dist(lambda x: pm.logsumexp(someDervDist))
p=0.05
p_vd = Mixture('p_vd',[p,1-p],[likelihoodDist1,likelihoodDist2],observed = actual_visit_mod_90_days)
错误信息是
---------------------------------------------------------------------------
IndexError Traceback (most recent call last)
<ipython-input-19-1c9d266976ce> in <module>()
28 p=0.05
29
---> 30 p_vd = Mixture('p_vd',[p,1-p],[likelihoodDist1,likelihoodDist2],observed = actual_visit_mod_90_days)
31
32 #print pm.distributions.generate_samples(postponeProbability)
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\distributions\distribution.pyc in __new__(cls, name, *args, **kwargs)
35 total_size = kwargs.pop('total_size', None)
36 dist = cls.dist(*args, **kwargs)
---> 37 return model.Var(name, dist, data, total_size)
38 else:
39 raise TypeError("Name needs to be a string but got: {}".format(name))
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\model.pyc in Var(self, name, dist, data, total_size)
780 var = ObservedRV(name=name, data=data,
781 distribution=dist,
--> 782 total_size=total_size, model=self)
783 self.observed_RVs.append(var)
784 if var.missing_values:
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\model.pyc in __init__(self, type, owner, index, name, data, distribution, total_size, model)
1228
1229 self.missing_values = data.missing_values
-> 1230 self.logp_elemwiset = distribution.logp(data)
1231 # The logp might need scaling in minibatches.
1232 # This is done in `Factor`.
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\distributions\mixture.pyc in logp(self, value)
111 w = self.w
112
--> 113 return bound(logsumexp(tt.log(w) + self._comp_logp(value), axis=-1).sum(),
114 w >= 0, w <= 1, tt.allclose(w.sum(axis=-1), 1),
115 broadcast_conditions=False)
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\pymc3\distributions\mixture.pyc in _comp_logp(self, value)
83 except AttributeError:
84 return tt.stack([comp_dist.logp(value) for comp_dist in comp_dists],
---> 85 axis=1)
86
87 def _comp_means(self):
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\theano\tensor\basic.pyc in stack(*tensors, **kwargs)
4583 dtype = scal.upcast(*[i.dtype for i in tensors])
4584 return theano.tensor.opt.MakeVector(dtype)(*tensors)
-> 4585 return join(axis, *[shape_padaxis(t, axis) for t in tensors])
4586
4587
C:\Users\mdevananda\Continuum\Anaconda2\lib\site-packages\theano\tensor\basic.pyc in shape_padaxis(t, axis)
4475 if not -ndim <= axis < ndim:
4476 msg = 'axis {0} is out of bounds [-{1}, {1})'.format(axis, ndim)
-> 4477 raise IndexError(msg)
4478 if axis < 0:
4479 axis += ndim
IndexError: axis 1 is out of bounds [-1, 1)
对我来说，问题似乎出在 DensityDist.dist 函数上。我也不确定我是否遗漏了什么。
编辑更新：考虑到事件 1 的历史模式 (actual_visit_mod_90_days)，我正在尝试了解事件 1 从开始日期起多少天后会发生。最初，我将 someDervDist 定义为一个等式：
someDervDist = dist2-dist1+dist3-1
给出了决定事件 1 的事件的各种概率之间的关系。
likelihoodDist2 = pm.DensityDist.dist(lambda x: pm.logsumexp(someDervDist))
会给出 logsumexp(someDervDist.logp()) 我没有写成 .logp 因为我认为它会给出对数概率（嗯，我在这里错了吗？）但是，我意识到这需要日志概率并将 someDervDist eqn 更改为 DensityDist 分布（我也尝试使用 Deterministic）。 然后，我假设这偏离事件 1 的周期性发生模型的概率很小 p=0.05。基本上我正在寻找类似 Mixture([p, 1-p], [ DiscreteUniform(lower = 0, upper = 89), someDervDist]).
更新的答案，毕竟使用 Constant 发行版不是一个好主意，它很快就会被弃用。所以也许一个小技巧会有所帮助，只需使用均值 someDervDist 和小方差的分布。
actual_visit_mod_90_days = [29, 76, 66, 66, 10, 63, 69, 75, 66, 87, 58]
with pm.Model() as pattern_study:
dist1 = pm.Binomial('dist1', n=14, p=0.5)
dist2 = pm.DiscreteUniform('dist2', lower = 0, upper = 89)
pp = pm.DensityDist('pp', jeffreys_logp, testval = 0.5)
dist3 = pm.Geometric('dist3', pp)
someDervDist = pm.Deterministic('someDervDist', dist2-dist1+dist3-1)
likelihoodDist2 = pm.DiscreteUniform.dist(someDervDist-1, someDervDist + 1)
likelihoodDist1 = pm.DiscreteUniform.dist(lower=0,upper=89)
p= 0.05
p_vd = pm.Mixture('p_vd', [p, 1-p],
[likelihoodDist1, likelihoodDist2],
observed = actual_visit_mod_90_days)
trace = pm.sample(1000, step=pm.Metropolis())
作为替代方案，您可以尝试使用模型的连续版本，如下所示：
with pm.Model() as pattern_study:
dist1 = pm.Normal('dist1', 7, 1.5)
dist2 = pm.Uniform('dist2', lower = 0, upper = 89)
dist3 = pm.Exponential('dist3', 1)
someDervDist = pm.Deterministic('someDervDist', dist2-dist1+dist3-1)
likelihoodDist2 = pm.Normal.dist(someDervDist, 0.1)
likelihoodDist1 = pm.Uniform.dist(lower=0,upper=89)
p = 0.05
p_vd = pm.Mixture('p_vd', [p, 1-p],
[likelihoodDist1, likelihoodDist2],
observed = actual_visit_mod_90_days)
trace = pm.sample(1000)